% function[variances] = compute_var(mlh,param,time,alt,range,type)
%This function computes the verical variance of the data for a given range.

% Inputs:
% - mlh: A vector of mixed layer heights to use as reference
% - var: A vector of data where it is necessary to compute the vertical
% variance
% -alt: A vector containing the altitudes
% - time: A vector of times for the data 
% - range: A range of altitude by which to compute the vertical variance
% (e.g., 200m,300m,400m)

% Outputs:
% - var: The vertical variance of the given altitude range


%--Get training and testing data from flight data--
trainFiles = {'20190708_R0.h5'	'20190726_R0.h5'	'20190627_R0.h5'	'20190725_R0.h5'	'20190709_R0.h5'	'20190720_R0.h5'	'20190727_R1.h5'	'20190617_R0.h5'	'20190710_R0.h5'	'20190622_R0.h5'	'20190707_R0.h5'};
testFiles = {'20190711_R0.h5'	'20190722_R0.h5'	'20190724_R0.h5'	'20190701_R0.h5'	'20190620_R0.h5'	'20190704_R0.h5'};

%--Extract full directories for training
arch_train = getFilePaths(trainFiles, '/Users/juliechristopoulos/Desktop/PBL/ACT_AMERICA/Archive');

%--Inputs--
[mlh] = gatherVar(arch_train, '/DataProducts/MixedLayerHeight');
[alt] = gatherVar(arch_train, '/DataProducts/Altitude');
[time] = gatherVar(arch_train, '/Nav_Data/gps_time');
[param] = gatherVar(arch_train,'/DataProducts/532_bsc');
range = 10;
type = 'horizontal';


if contains(type,'vertical')
    variances = cell(size(mlh));
    for vector = 1:length(mlh)
        mlh_data = mlh{vector};
        mlh_data(isnan(mlh_data)) = 0;
        alt_data = alt{vector};
        time_data = time{vector};
        param_data = param{vector};
        variances{vector} = [];
        for i = 1:length(time_data)
            mlh_time = mlh_data(i);
            param_time = param_data(:,i);
            if (mlh_time ~=0)
                min_index = fix(interp1(alt_data,find(alt_data),mlh_time-range));
                max_index = fix(interp1(alt_data,find(alt_data),mlh_time+range));
                param_range = param_time(min_index:max_index);
                param_variance = var(param_range,"omitnan");
                variances{vector} = [variances{vector}, param_variance]; % store value into cell array
            else
                variances{vector} = [variances{vector}, nan]; % store NaN value into cell array
            end
        end
    end
end 

if contains(type,'horizontal')
    variances = cell(size(mlh));
    for vector = 1:length(mlh)
        mlh_data = mlh{vector};
        mlh_data(isnan(mlh_data)) = 0;
        alt_data = alt{vector};
        time_data = time{vector};
        param_data = param{vector};
        variances{vector} = [];
        for i = 1:length(time_data)
            mlh_time = mlh_data(i);
            if (mlh_time ~= 0)
                alt_index = fix(interp1(alt_data,find(alt_data),mlh_time));
                i_start = i - range;
                i_end = i + range;
                param_range = param_data(alt_index,i_start:i_end);
                param_variance = var(param_range,"omitnan");
                variances{vector} = [variances{vector}, param_variance]; 
                
            else
                variances{vector} = [variances{vector}, nan]; 
            end
        end
    end

end



    
    













