% CPC TSI data processing - T1 UEA
% Written by Luciana Rizzo
% Flex data files (Labview DAQ)
% Read raw data files, apply corrections, calculate averages, make plots
% Functions used: frepeated; export_LFA_GoAmazon_jan2015
% Outputs: cpc_Level0 (mattime(utc);conc(cm-3))--> all raw data, without filters
%          cpc_Level1_5min (dd-mmm-yyyy HH:MM; sec_since_1970; mean; std; median; p25; p75; max; min; num)
%          cpc_Level1_30min
% Time stamps in the output files refer to the start of the period  
%
% Level 1: 1-min conc must be in the range 50-10^5 cm-3. No instrument error should be reported.
% Level 2: Apply STP corrections and coincidence corrections.
%          Os dados flex vêm direto da saída serial do CPC, utilizando o comando RALL.
%          Do manual do CPC, entendi que a concentração que sai via serial já está corrigida para erros de coincidência. 
% Level 3: Intercomparison with other CPC if available. Particle counts comparison must be within 20%.
%
% Parametros de entrada --------------------------------------------
caminho1='C:\A\GOAMAZON\T1\CPC_3776151301sn\csv\'; % arquivos flex (Labview)
caminho2='C:\A\GOAMAZON\T1\CPC_3776151301sn\';
tit='T1_2015-2016';
pref=1013.25; % Reference pressure for STP correction (mbar)
tref=293.15; % Reference temp for STP correction (K)
corrigirtempo=false;
filtmanut=false; % logical para filtrar dados em períodos de manutenção ou com problemas
readprev=false; % logical for reading previous processed data 
lev1=false;
lev2=false;
%
% inicio e termino de periodos de manutencao (horario em utc): removidos das analises
% if(filtmanut)
%     C1=xlsread('C:\A\Diadema2016\logbook_trailer_Diadema2016.xlsx','todos','A4:M135');
%     C2=xlsread('C:\A\Diadema2016\logbook_trailer_Diadema2016.xlsx','cpc','A4:M13');
%     manut(:,1)=cat(1,datenum(C1(:,1:6)),datenum(C2(:,1:6)));
%     manut(:,2)=cat(1,datenum(C1(:,7:12)),datenum(C2(:,7:12)));
%     manut(:,3)=cat(1,C1(:,13),C2(:,13));    
% end
%
% Program start -----------------------------------------------------
files=dir(caminho1);
alldata=zeros(1,2);
if (readprev)
    prevdata=load([caminho2,tit,'_cpc3772_level0.mat']);
    lerdaqui=floor(prevdata.alldata(length(prevdata.alldata),1));    
else
    lerdaqui=datenum(1900,1,1);
end
disp('-----Reading flexdata files')
for ifile=3:size(files,1)
    nome=getfield(files, {ifile,1}, 'name');
    filedia=datenum(nome(length(nome)-13:length(nome)-4),'yyyy-mm-dd');
    if(filedia>=lerdaqui)
        disp(nome)
        fid=fopen(strcat(caminho1,nome));
        C=cell(floor(getfield(files, {ifile,1}, 'bytes')/92),3);
        i=1;
        while (~feof(fid))
            C(i,:)=textscan(fid,'%*s%s%s %f %*s%*s%*s%*s%*s%*s%*s%*s%*s%*s%*s',1,'Delimiter',',');
            i=i+1;
        end
        fclose(fid);
        %
        % removendo linhas com defeito ou cabeçalhos
        rem=0;
        for j=3:3
            aux=find(cellfun(@isempty,C(:,j)));
            rem=cat(1,rem,aux);
        end
        rem(1)=[];
        C(rem,:)=[];
        rem=find(cellfun(@isnan,C(:,3)));
        C(rem,:)=[];
        %
        data_(:,2)=cell2mat(C(:,3));
        for i=1:length(C)
            data_(i,1)=datenum(strcat(C{i,1},C{i,2}),'yyyy-mm-ddHH:MM:SS');
        end
        alldata=cat(1,alldata,data_);
        clear data_ C
    end    
end
alldata(1,:)=[];
if (readprev)
    aux=find(prevdata.alldata(:,1)<min(alldata(:,1)));
    alldata=cat(1,prevdata.alldata(aux,:),alldata);
end
alldata=sortrows(alldata,1); % sort data based on day of year
%
if(corrigirtempo) % correção segundo o logbook
end
% 
% Obs: é comum a repetição de horários nos arquivos flexdata. Nada grave,
% dura alguns segundos por dia, com valores de concentração ligeiramente diferentes.
% Tem a ver com um descompasso entre o instante de coleta do LabView e a
% frequência de medidas do CPC. Assim, vou dispensar a verificação de 
% dados repetidos aqui, e simplesmente remover linhas com mesmo time stamp.
disp('----Checking for repeated data')
aux=find(diff(alldata(:,1))==0);
disp(['Número de linhas com timestamp repetido: ',num2str(length(aux))])
alldata(aux,:)=[];
clear aux
% [newdata,erro]=frepeated(alldata,2);
% if (strcmp(erro,'none')~=1) 
%     disp(erro)
% else
%     alldata=newdata;
% end
% clear aux erro 
%
save([caminho2,tit,'_cpc3776_level0.mat'],'alldata','-mat') 
% 
% Remove períodos de manutencao/contaminacao de acordo com o logbook
% if(filtmanut)
%    alldataorig=alldata;
%    for i=1:size(manut,1)
%         aux=find(alldata(:,1)>=manut(i,1) & alldata(:,1)<=manut(i,2));
%         if(manut(i,3)==0)            
%             alldata(aux,:)=[]; % remove maintanance periods from dataset
%         end    
%    end
%    save([caminho2,tit,'_cpc3772_clean.mat'],'alldata','-mat') 
%    csvwrite([caminho2,tit,'_cpc3772_clean.csv'],[alldata(:,1:2),datevec(alldata(:,1))]) 
%    figure,plot(alldataorig(:,1),alldataorig(:,2),'r.'),hold on,grid on
%    plot(alldata(:,1),alldata(:,2),'c.')
%    ylabel('CPC3772 (cm^-^3)'),title(tit)
%    dynamicDateTicks   
%    saveas(gcf,strcat(caminho2,tit,'_cpc3772_level0'),'jpg')
% else   
%     figure
%     plot(alldata(:,1),alldata(:,2),'-r.'),grid on,hold on
%     ylabel('CPC3776 (cm^-^3)'),title(tit)
%     set(gca,'XLim',[datenum(2015,12,2) datenum(2016,4,26)])
%     dynamicDateTicks
% end
%
%
% ------------ Level 1 ------------------------------------------------
% Remove local contamination (clean data)
%
if(lev1)
    %
    origsize=size(alldata,1);
    orig=alldata;
    rem=find(alldata(:,2)<50 | alldata(:,2)>10^5);
    alldata(rem,:)=[]; 
    %
    figure,plot(orig(:,1),orig(:,2),'-r.'),hold on,grid on
    plot(alldata(:,1),alldata(:,2),'b.'),ylabel('CPC3776 (cm^-^3)'),title(tit),dynamicDateTicks
    %
    % Save Level 1 data
    Instrument='CPC';
    Model='3776';
    serial_num='3776151301';
    Site=tit;
    Val_level='1';
    Proc_vers='1.0';
    labels{1}='conc';
    units{1}='(cm-3)';
    export_LFA_GoAmazon_jan2015(labels,units,Instrument,Model,serial_num,Site,Val_level,Proc_vers,origsize,alldata(:,1),alldata(:,2))
end
%
% ------------ Level 2 ------------------------------------------------
% Level 2: Apply STP corrections and coincidence corrections.
% Os dados flex vêm direto da saída serial do CPC, utilizando o comando RALL.
% Do manual do CPC, entendi que a concentração que sai via serial já está
% corrigida para erros de coincidência. 
if(lev2)        
    % Lê os arquivos Level 1 
    aux=importdata('C:\A\GOAMAZON\T1\LFA_standard\T1_2015-2016_CPC_3776_Level1_5min.csv');
    level1_5min=aux.data(:,1:2);        
    level1_5min(:,1)=level1_5min(:,1)/(24*60*60)+datenum(1970,1,1); % unix time to matlab time
    origsize=sum(~isnan(level1_5min(:,2)));
    orig=level1_5min;
    tit='T1_2015-2016';
    %
    % Lê fator stp
    C=importdata('C:\A\GOAMAZON\T1\LFA_standard\T1_2015-2016_AWS_Gill_Level2_5min.csv');
    stp=[C.data(:,1) C.data(:,10)];  
    stp(:,1)=stp(:,1)/(24*60*60)+datenum(1970,1,1); 
    level1_5min(:,2)=level1_5min(:,2).*stp(:,2);
    level1_5min(level1_5min(:,1)>datenum(2016,3,13),:)=[];
    %
    figure, plot(orig(:,1),orig(:,2),'r.'),grid on, hold on,dynamicDateTicks
    plot(level1_5min(:,1),level1_5min(:,2),'b.')
    %
    % Save Level 2 data 
    Instrument='CPC';
    Model='3776';
    serial_num='3776151301';
    Site=tit;
    Val_level='2';
    Proc_vers='1.0';
    labels{1}='conc';
    units{1}='(cm-3)';
    export_LFA_GoAmazon_jan2015(labels,units,Instrument,Model,serial_num,Site,Val_level,Proc_vers,origsize,level1_5min(:,1),level1_5min(:,2))
end
%
%
% ------------ Level 3 ------------------------------------------------
% Level 3: Compare to SMPS integrated concentration 
if(lev3)        
    % Lê o arquivo Level 2 
    aux=importdata('C:\A\GOAMAZON\T1\LFA_standard\T1_2015-2016_CPC_3776_Level2_5min.csv');
    level2=aux.data;
    level2(:,1)=level2(:,1)/(24*60*60)+datenum(1970,1,1); % unix time to matlab time    origsize=sum(~isnan(level1_5min(:,2)));
    orig=level2;
    origsize=sum(~isnan(level2(:,2)));
    tit='T1_2015-2016';
    %
    aux=importdata('C:\A\GOAMAZON\T1\LFA_standard\T1_2015-2016_SMPS_TSI_3081_Level3_5min.csv');
    param=str2double(aux.textdata(3:length(aux.textdata),2:5)); % TotConc(cm-3); TotVol(cm3); Dpg(nm)
    param(:,1)=param(:,1)/(24*60*60)+datenum(1970,1,1); % unix time to matlab time
    razao=param(:,2)./level2(:,2);
    %
    figure
    subplot(2,1,1),plot(param(:,1),param(:,2),'b.'),grid on,hold on,dynamicDateTicks
    plot(level2(:,1),level2(:,2),'r.')
    ylabel('cm^-^3','FontSize',12),title(tit,'FontSize',12)
    subplot(2,1,2),plot(param(:,1),razao,'b.'),grid on,dynamicDateTicks
    ylabel('Ratio SMPS/CPC','FontSize',12),title(tit,'FontSize',12)
    %
    level2(razao>=2.5,2)=NaN;
    level2(level2(:,2)<=200,2)=NaN;
    %
    % Save Level 2 data 
    Instrument='CPC';
    Model='3776';
    serial_num='3776151301';
    Site=tit;
    Val_level='3';
    Proc_vers='1.0';
    labels{1}='conc';
    units{1}='(cm-3)';
    export_LFA_GoAmazon_jan2015(labels,units,Instrument,Model,serial_num,Site,Val_level,Proc_vers,origsize,level2(:,1),level2(:,2))    
end






