function [H,Hn,ny,Tc] = RAPTOR_diag_sim(X,G,V,U,model)
% simulate diagnostics output based on model.diagnistics.diagfun
% function [H,Hn,ny,Tc] = RAPTOR_diag_sim(X,G,V,U,model)
% H: clean measurements
% Hn: measurements corrupted by noise
% ny: number of diagnostic channels
% scaling matrix from H to physical measurements

assert(isa(model.diagnostics.diagfun,'function_handle'),'model.diagnostics.diagfun has to be a function handle')

diagfun = model.diagnostics.diagfun;
for it = 1:size(X,2)
    [hk,~,~,Rk,Tc,iTc] = diagfun(X(:,it),G(:,it),V(:,it),U(:,it),it,model);
    v = sqrt(Rk)*randn(size(hk)); % noise
    H(:,it) = hk;
    Hn(:,it) = H(:,it) + v;
end
ny = size(H,1);
