%% Read lamp 

addpath(genpath('C:\Users\ecimoli\OneDrive\Documenti\MATLAB')); 
addpath(genpath('E:\2018_Cape_Evans_data'));

clear
close all
cd 'E:\2019_Cameras_calibration\ResononPikaII_OCT2019\RAD_cal_Oct_2019\f2_0focus0_46m';

% Choose lamp file & folder 
% Kestrel f2.8
% lampD_Kestrel = readtable('lampD_attenuator100open_AisaKestrel10.csv');
% 
% RAD_lamp = lampD_Kestrel.resampledSpectralValue_W_m2_sr_nm_;

% Pika IIi lampA f1.4
% 
% lampA_Pika = readtable('lampA_attenuator100open_ResononPika2.csv');
% 
% RAD_lamp = lampA_Pika.resampledSpectralValue_W_m2_sr_nm_;
% 
% % Pika IIi lampD f2.0
% 
lampD_Pika = readtable('lampD_attenuator100open_ResononPika2.csv');

RAD_lamp = lampD_Pika.resampledSpectralValue_W_m2_sr_nm_;

% Read RAD & DC

info=envihdrread('RAD_IT_05ms.hdr');
[HI_RAD]=envidataread('RAD_IT_05ms.bil',info);
[x_RAD,y_RAD,bands_RAD] = deal(info.lines, info.samples, info.bands);

info=envihdrread('DC_f2_IT9ms.hdr');
[HI_DC]=envidataread('DC_f2_IT9ms.bil',info);
[x_DC,y_DC,bands_DC] = deal(info.lines, info.samples, info.bands);

info.wavelength(1) = ''; % change caps for kestrel (W)/pika (w)
info.wavelength = info.wavelength(1:end-1);
WVrange = strsplit(info.wavelength,',');
WVrange = str2double(WVrange(:,:)');

% IT = info.tint; %for kestrel
IT = info.shutter; %for pika

% Average along track

mean_HI_RAD = squeeze(mean(HI_RAD,1));

mean_HI_DC = squeeze(mean(HI_DC,1));

% Compute coefficient

RAD_lamp_matrix = repmat(RAD_lamp,1,info.samples);
RAD_lamp_matrix = RAD_lamp_matrix';

coeff_HI = RAD_lamp_matrix./(mean_HI_RAD - mean_HI_DC) * IT;
%% open folder
cd 'E:\2019_Cameras_calibration\ResononPikaII_OCT2019\Coeff_Pika_f2_0_2x1';

%% Convert into seconds

coeff_HI = coeff_HI/1000;

plot(WVrange,coeff_HI','DisplayName','coeff_HI')

%% Write 

cd 'E:\2019_Cameras_calibration\ResononPikaII_OCT2019\Coeff_Pika_f2_0_2x1_sec';

save('PikaII_IT9_2_0','coeff_HI');

csvwrite('PikaII_IT9_2_0.csv',coeff_HI);

%% Plot

plot(coeff_HI','DisplayName','coeff_HI');

