%--------------------------------------------------------------------
%Normalization and correction for water clusters according to 
% Warneke & DeGouw (2007) - Mass Spec Rev, 26, 223-257
fid=fopen('PTRMS/waterco.txt');
waterco=textscan(fid,'%d %d','Delimiter',' ');
fclose(fid);
waterco=cell2mat(waterco); %cell2mat - concatenation of ??
waterco=double(waterco);

%waterco.txt has Xr * 100 (as proposed by Warneke & DeGouw) that's why here
%it is multiplied by 0.01:
waterco(:,2)=0.01.*waterco(:,2);

dummy=zeros(1,max(size(header))); %creates an array of zeros of 1 row and n of header columns

% Alimenta os mz corretos com os valores de waterco que existem 
for i=1:max(size(waterco))
    [diff,idx]=min(abs(waterco(i,1)-header));
    %diff -> o valor da operacao 'min(abs(waterco(i,1)-header))'
    %idx -> a posicao no vetor em que diff foi encontrado
    %quando diff = 0, significa que waterco contem o valor de mz presente
    %no header. Entao, dummy(idx) recebe o valor da correcao (colum 2 of
    %waterco).
    if diff==0
        dummy(idx)=waterco(i,2);
    end;
end

%atulaiza waterco
waterco=dummy;

%--------------------------------------------------------------------
%Here it starts the calculation of mixing ratios 

mat_Prim=Prim*ones(1,max(size(header)));  

mat_Prim=mat_Prim';

% correction for water clusters
mat_m37=waterco'*m37';

mat_waterco=1e6./(mat_Prim+mat_m37)';

%'vocs' are the mixing ratios for each VOC
vocs=vocs.*mat_waterco;

%--------------------------------------------------------------------

%First IOP (intensive operation period)        
factor_31=1; 
factor_33=14.69;
factor_42=25.93;
factor_45=19.93;        
factor_47=3; %Factory calibration factor table
factor_59=26.2;
factor_61=5.6; 
factor_63=18.16;
factor_69=12.02; 
factor_71=14.4; 
factor_73=24.84; 
factor_79=17.14; 
factor_93=19.42; 
factor_107=19.62; 
factor_121=18.04;
factor_137=7.09; 

% check if there is a Cal_Factor for the respective mz. 
for i=1:size(header,2)
    eval(['dummy=exist(''factor_',num2str(header(i)),''',''var'');'])
    %If not, use Factor = 10
    %@Milena: why 10??
    if dummy==0
        eval(['factor_',num2str(header(i)),'=10;'])
    end;
    eval(['cal_factor(i)=factor_',num2str(header(i)),';'])
end;

cal_factor=ones(max(size(time_PTRMS)),1)*cal_factor;

%Second IOP - 560V DT        
% factor_31=1.19; 
% factor_33=13.409;
% factor_42=25.26;
% factor_45=15.668;        
% factor_47=3; %Factory calibration factor table
% factor_59=30.06;
% factor_61=5.6; 
% factor_63=17.19;
% factor_69=13.94; 
% factor_71=25.74; 
% factor_73=27.42; 
% factor_79=16.45; 
% factor_93=18.71; 
% factor_107=19.47; 
% factor_121=13.72;
% factor_137=8.01; 
         


for i=1:size(header,2)
    eval(['cal_factor(time_PTRMS>datenum(2014,07,30),i)=factor_',num2str(header(i)),';'])
end;



%Campaign T1
factor_31=0.50; 
factor_33=9.86;
factor_42=18.64;
factor_45=15.93;        
factor_47=3; %Factory calibration factor table
factor_59=30.06;
factor_61=5.6; 
factor_63=21.63;
factor_69=8.658; 
factor_71=17.47; 
factor_73=21.5349; 
factor_79=11.65; 
factor_93=20.577; 
factor_107=14.5094; 
factor_121=13.72;
factor_137=5.00102; 
         

for i=1:size(header,2)
    eval(['cal_factor(time_PTRMS>datenum(2015,12,01),i)=factor_',num2str(header(i)),';'])
end;


%Campaign ZF2
factor_31=6.92; 
factor_33=4.44;
factor_42=17.84;
factor_45=14.15;        
factor_47=3; %Factory calibration factor table
factor_59=18.11;
factor_61=5.6; 
factor_63=21.63;
factor_69=11.869; 
factor_71=30.53; 
factor_73=12.6133; 
factor_79=11.59; 
factor_93=34.51; 
factor_107=14.5094; 
factor_121=7.57;
factor_137=3.981; 
         
for i=1:size(header,2)
    eval(['cal_factor(time_PTRMS>datenum(2016,06,01),i)=factor_',num2str(header(i)),';'])
end;


vocs=vocs./cal_factor;


%label each mz
for i=1:size(header,2)
    eval(['label_',num2str(header(i)),'=[label_',num2str(header(i)),', '' (ppbv)''];']);
end;


