function lic = read_licor_CAPRICORN_2016(dfl,ddd,hhh)
%{
reads son file specified by dfl and returns array with the
following columns:

1   jd_ref      10 Hz timestamp
2   Licor_H2O   H2O, mmol/m3
3   Licor_CO2   CO2, approx mmol/m3
4   Licor_CO2a  CO2, approx mol/mol
5   Licor_T     Licor box temperature
6   Licor_P     Licor box pressure kPa
7   AGC         Licor diagnostic

input parameters: dfl = string path to scs file
                  ddd = julian date
                  hhh = hour
%}
%% initialize
disp(['Reading lic file for hour ',int2str(hhh)]);
flist = fopen(dfl,'r');
temp = {};
z = [];

%% read data
while feof(flist)==0,
    try    % read file into cell array. temp is 1x10 cell array
        temp = textscan(flist,'%2f%2f%3f %f %f %f %f %f %f %f %*[^\n]','delimiter', ', ','headerlines', 1,'emptyvalue',NaN,'treatAsEmpty','Overflow');
        z = [z cell2mat(temp)'];
    catch
        for ii=1:10  % length of last cell reflects missing values from any/all fields
            if length(temp{1,ii})~=length(temp{1,10})
                temp{1,ii}(length(temp{1,10})+1)=[]; %truncate length of all cells to length(temp{1,10})
            end;
        end;
        if ~isempty(temp{1,1})
            z=[z cell2mat(temp)'];
        end;
    end;
end;

z(:,any(isnan(z),1))=[]; % kill lines beginning with NaN

chk = size(z);
if chk(1)==10    % licor set up to output 10 values...
    Licor_diag = z(4,:);
    Licor_CO2_Abs = z(5,:); % CO2 Absorbance
    Licor_CO2 = z(6,:); % mmol/m^3
    Licor_H2O_Abs = z(7,:); % H2O Absorbance
    Licor_H2O = z(8,:); % mmol/m^3
    Licor_T = z(9,:); % C
    Licor_P = z(10,:); % kPa
end;
fclose(flist);
tlic =  ddd + (hhh + (z(1,:)+(z(2,:)+z(3,:)/1000)/60)/60)/24;

%% convert CO2 units from mol/m3 to best estimate of mol/mol
Rgas_universal = 8.314472; % Pa m3 /K /mol
Licor_CO2a = (Licor_CO2*1e-3).*Rgas_universal.*(Licor_T+273.15)./(Licor_P*1e3); % mmol/m3 to mol/mol

%% interpret status diagnostic
xx = Licor_diag;
b7 = floor(xx/128);
y = b7*128;
b6 = floor((xx-y)/64);
y = y+b6*64;
b5 = floor((xx-y)/32);
y = y+b5*32;
b4 = floor((xx-y)/16);
y = y+b4*16;
AGC = (xx-y)*6.25;

%% interpolate to 10 Hz reference timestamp
% No output for absorbances in this version...
start = (ddd+hhh/24);
delta = double(1.0/864000);
last = start + 35999*delta;
tref = start:delta:last;
lic = zeros(36000,7);
lic(:,1) = tref;
lic(:,2) = interp1(tlic,Licor_H2O,tref,'linear','extrap');
lic(:,3) = interp1(tlic,Licor_CO2,tref,'linear','extrap');
lic(:,4) = interp1(tlic,Licor_CO2a,tref,'linear','extrap');
lic(:,5) = interp1(tlic,Licor_T,tref,'linear','extrap');
lic(:,6) = interp1(tlic,Licor_P,tref,'linear','extrap');
lic(:,7) = interp1(tlic,AGC,tref,'linear','extrap');


