% routine to 
% 1.) load raw data
% 2.) create a time grid with 0.1 s spacing
% 3.) flip sign of u and v
% 4.) rotate data into meteorological coordinates
% 5.) store into matlab mat file

clear all;
close all;

searchpath='./raw/';
target_directory='./mat/';
filenameextension='*.dat';
files=dir(strcat(searchpath,filenameextension));

for i=1:length(files)
    filename=files(i,1).name % run any code you want to have applied to all files

    %**************************************************************************
    % read in the fast way
    try
        [yyyy,mm,dd,hh,minmin,sec,u,v,w,temp,errorcode,battery] = ...
        textread(strcat(searchpath,...
                        filename),...
                        '%u %u %u %u %u %f   %f %f %f %f  %2c %f');
        serial_date_number = datenum(yyyy,mm,dd,hh,minmin,sec);
        errorcode=hex2dec(char(errorcode));
        clear('yyyy','mm','dd','hh','minmin','sec');
    catch
    % read in the slow way if it fails
    % retrieve number of data lines to be read:
    tic;
        N = linecount(strcat(searchpath,filename));
    toc; 
    % Open ASCII file:
    fid=fopen(strcat(searchpath,filename));
                
    % Pre-allocate arrays for wind velocity components (u, v, w), 
    % temperature in deg Celsius (Tc), and MATLAB serial time (t).
    u = zeros(N,1);
    v = zeros(N,1);
    w = zeros(N,1);
    temp = zeros(N,1);
    serial_date_number = zeros(N,1);
    errorcode = zeros(N,1);
    battery = zeros(N,1);

    % Read header (first two lines):
    %tline = fgetl(fid)
    %tline = fgetl(fid)

    % Read data line by line as text 
    % and convert to double-precision numbers:
    
    tic;
    for i = 1:N
        tline   = fgetl(fid);
        year    = str2num(tline(1:4));
        month   = str2num(tline(6:7));
        day     = str2num(tline(9:10));
        hour    = str2num(tline(12:13));
        minute  = str2num(tline(15:16));
        second  = str2num(tline(18:26));
        u(i)    = str2num(tline(27:34));
        v(i)    = str2num(tline(35:41));
        w(i)    = str2num(tline(42:48));
        temp(i)   = str2num(tline(49:56));
        errorcode(i)   = hex2dec(tline(57:60));
        battery(i)   = str2num(tline(61:end));
        serial_date_number(i)    = datenum(year,month,day,hour,minute,second);
    end
    toc;
    end
%**************************************************************************

%**************************************************************************
% save raw - not averaged data
mat_filename = strcat(target_directory,...
                      filename(1:8),...
                      'u_v_w_temp_',...
                      datestr(serial_date_number(1),30),...
                      '_',...
                      datestr(serial_date_number(end),30),...
                      '.mat');
tic;             
    save(mat_filename,'serial_date_number','u','v','w','errorcode','temp','battery','-mat');   
toc;
%**************************************************************************

%**************************************************************************
% create grid for given time interval
grid_interval=0.1 % sec;
begin_serial_date_number=serial_date_number(1);
end_serial_date_number=serial_date_number(end);
number_of_avg_intervals=fix((end_serial_date_number-begin_serial_date_number).*(24*3600)./grid_interval);

% first grid element
first_grid_date_num=datevec(serial_date_number(1));
tens=mod(first_grid_date_num(6),fix(first_grid_date_num(6)))./grid_interval;
new_seconds_grid=fix(first_grid_date_num(6))+(fix(tens).*grid_interval);
first_grid_serial_date_number=datenum(first_grid_date_num(1),...
                                      first_grid_date_num(2),...
                                      first_grid_date_num(3),...
                                      first_grid_date_num(4),...
                                      first_grid_date_num(5),...
                                      new_seconds_grid);
                                  
% print the first gridded serial date number in a human readable format
% datevec(first_grid_serial_date_number)
% create rest of the grid elements
grid_serial_date_number=first_grid_serial_date_number+datenum(0,0,0,0,0,[0:number_of_avg_intervals].*grid_interval);
%**************************************************************************
% preallocate array for gridded data
u_gridded = zeros(number_of_avg_intervals,1);
v_gridded = zeros(number_of_avg_intervals,1);
w_gridded = zeros(number_of_avg_intervals,1);
temp_gridded = zeros(number_of_avg_intervals,1);
serial_date_number_gridded = grid_serial_date_number';
errorcode_gridded = zeros(number_of_avg_intervals,1);
battery_gridded = zeros(number_of_avg_intervals,1);
%**************************************************************************
% BEGIN actuall averaging into grid structure
%**************************************************************************
tic;
for j=1:(number_of_avg_intervals)
    %  j
    %  datevec(serial_date_number_gridded(j))
    
    upper_bound=(serial_date_number_gridded(j)+datenum(0,0,0,0,0,grid_interval/2));
    % datevec(upper_bound)
    
    lower_bound=(serial_date_number_gridded(j)-datenum(0,0,0,0,0,grid_interval/2));
    % datevec(lower_bound)
    
    % serial_date_matching_index_lower=find(serial_date_number > lower_bound);
    % serial_date_matching_index_upper=find(serial_date_number < upper_bound);

    % to find matching serial date number, that is at least on original data
    % point is available within the grid interval
    try
        serial_date_matching_index=find(((serial_date_number > lower_bound)&(serial_date_number < upper_bound)));
        u_gridded(j)=mean(u(serial_date_matching_index));
        v_gridded(j)=mean(v(serial_date_matching_index));
        w_gridded(j)=mean(w(serial_date_matching_index));
        temp_gridded(j)=mean(temp(serial_date_matching_index));
        battery_gridded(j)=mean(battery(serial_date_matching_index));
        errorcode_gridded(j)=max(errorcode(serial_date_matching_index));
    catch
    % go here if no data available for a given grid element
        u_gridded(j)=NaN;
        v_gridded(j)=NaN;
        w_gridded(j)=NaN;
        temp_gridded(j)=NaN;
        battery_gridded(j)=NaN;
        errorcode_gridded(j)=NaN;
    end
end
toc;
%**************************************************************************
% END actuall averaging into grid structure
%**************************************************************************

%**************************************************************************
% free up some memory
clear('mat_filename','serial_date_number','u','v','w','temp','battery','errorcode');
%**************************************************************************

mat_filename_gridded = strcat(target_directory,...
                      filename(1:8),...
                      'u_v_w_temp_gridded_',...
                      datestr(serial_date_number_gridded(1),30),...
                      '_',...
                      datestr(serial_date_number_gridded(end),30),...
                      '.mat');
tic;             
    save(mat_filename_gridded,'serial_date_number_gridded','u_gridded','v_gridded','w_gridded','errorcode_gridded','temp_gridded','battery_gridded','-mat');   
toc;

%**************************************************************************
% free up some memory
clear('mat_filename_gridded','serial_date_number_gridded','u_gridded','v_gridded','w_gridded','temp_gridded','battery_gridded','errorcode_gridded');
%**************************************************************************
end