%CONVERTLIDAR  script in progress
%
%See also: 

% function convertLidar()
% Insert header doc here....

% This script can be used to convert lidar transect (raai) files into
% netcdf

% A few entities are used here:
% lidar transects: as read from the lidar transect file
% transect: assumed to be measurements along a seaward angle
% grid: in this case a curvilinear structured grid 
% netcdf output: a netcdf file following the cf convention where results
% are stored.
% 
% Transect
% A transect is considered to be a set of measurements along a line. This
% line has a local seaward coordinate 

% Lidar transect
% USGS performs lidar measurements resulting in a cloud of points along
% US coast. Each datapoint is assigned to a transect defined along the
% coast. The data is further devided in regions in the US.
% Each raai has a region specific distance within the region along the
% coastline.

% Grid
% The grid which contains the data for the transects is curvilinear
% structured. Each point in the grid has a seperate geographic coordinate.
% The grid is considered to be curvilinear, it follows the coastline. 
% The grid is structured, all points are related to 2 dimensions (seaward,
% along the coast)

% Netcdf output
% The netcdf file is stored according to the CF convention

% This code was written by Ben deSonnivlle in July of 2010
% afarris@usgs.gov added many comments in May of 2010
%       
% Let me (afarris) add a bit to Ben's comments above 
%  What Benn calls transects are our cut profiles.  They have many data
%  points per profile.  I think of the primary data being shoreline
%  position which varies as a function of distance along reference line (and
%  time (you can think of the data as a 3D matrix).  This organization
%  makes it difficult to store the profile data with the shoreline 
%  (although it can be done in Matlab with a structure).  Ben solved this
%  problem by creating a grid.  This seems to solve the storage problem
%  nicely, but I don't normally think of the data this way so I found it a
%  bit confusing.  Ben refers to the shoreline and bias data as meta-data,
%  altho' I don't think he stores it as metadata in the netcdf file (but maybe
%  I am wrong).  I don't agree with this view, to me it is the primary
%  data....but I can see how one could view these single data points per
%  profile as being meta-data to the main profiles (many data points per
%  profile). 

% the terms transect and profile are used interchangably

% Input:
% Data is stored in https:

%% Rekenpc
%  addpath('P:\mctools\matlab\io\nctools\mexnc\');
%  addpath('P:\mctools\matlab\io\snctools\');
 
%% Define path settings

dataType = 'Transects';
dataSet = 'Lidar_USGS';
rawDataDir = '..\raw\';
metaDataDir = '..\raw\';
nc_filename = 'washington_test.nc';
rawFileArray = dir([rawDataDir, '*.txt']);   
dbstop if error
tic
%% Allocate structure for dataset. 
% I think this statement is reduntant, since 'updateLidarGrid.m' calls it too.  
% However, since the code runs, I am not messing with it. (af2010May17)
transectStruct = createlidartransectstruct();
transectStruct(:) = [];

%% Load transect (profile) data
% Loop over all *.txt files
% Concatenate them into a big transectStructstructure Array.
% This loads the cut profile data (many data points per profile).
% It creates the variables:
% transect.lon                      longitude of each lidar data point
% transect.lat                      latiitude of each lidar data point
% transect.alititude                height of each lidar data point
% transect.crossShoreCoordinate     dist. from data point to profile origin
% transect.id                       profile number
% transect.easting                  easting of each lidar data point
% transect.northing                 northing of each lidar data point
for i = 1:length(rawFileArray)
    filename = [rawDataDir filesep rawFileArray(i).name];
    disp(['reading ' filename])   
    transectStruct = [transectStruct lidar2transect(filename, rawFileArray)];
end

% define a grid based on the transectStruct
% the data is transformed from things like transect.lon to grid.lon
% (other things are done too)
grid = transect2lidargrid(transectStruct);

%% Load shoreline and bias data, add to profile data
% (from *.mat file) 
% Ben often calls this the meta data.  There is one data point per profile.
% These are the variables that are added:
% grid.id                       profile number
% grid.shorePos                 distance from shoreline point to profile origin
% grid.shoreNorth               northing of shoreline position
% grid.shoreEast                easting of shoreline position
% grid.actionCode               this code tells us more about data 
% grid.bias                     
% grid.significant_wave_height 
% grid.deep_water_wave_length 
% grid.beach_slope 

% TODO add more variables...

grid = updateLidarGrid(grid, metaDataDir); 

%% store grid and transects to a netcdf file
% this code actually writes the netcdf file
lidargrid2netcdf(nc_filename, grid);

%% Determine whether all data can be written in one go

if sum([rawFileArray.bytes]) < 1E9
    lidartransect2netcdfblock(nc_filename, transectStruct);
else
    for j = 1:length(rawFileArray)
        filename = [rawDataDir filesep rawFileArray(j).name];
        disp(['reading data of ' filename]);
        transectStruct = [lidar2transect(filename, rawFileArray(j))];
        lidartransect2netcdf(nc_filename, transectStruct);
    end
end
    
disp('Done converting lidar data')
toc

% Elapsed time is 31684.693122 seconds.

