source: trunk/src/read_lai.m @ 327

function [lon_value, lat_value, datestr_value, lai_value]=read_lai(yyyy)

%READ_LAI fill array LAI 3d array and associated arrays date, lat, lon
%          with data in :file:`${PROJECT_ID}/LAI/laisen{yyyy}_float.txt`

%
%+
%
% ==========
% read_lai.m
% ==========
%
% .. function:: read_lai(yyyy)
%
% DESCRIPTION
% ===========
%
%    :param yyyy: year
%    :type yyyy: int16
%    :raise yyyy: required
%
% read a file :file:`${PROJECT_ID}/LAI/laisen{yyyy}_float.txt`
%
% Les LAI sont stockées dans le tableau 3D (temps, lat, lon) **lai_value**
%
% Les longitudes sont stockées dans le tableau 1D **lon_value**.
%
% Les latitudes sont stockées dans le tableau 1D **lat_value**.
%
% Les dates au format `yyyymmdd` sont stockées dans le tableau 1D **datestr_value**.
%
% EXAMPLES
% ========
%
% To read :file:`${PROJECT_ID}/LAI/laisen2001_float.txt`::
%
%   $ octave
%   octave> varamma_startup
%   octave> more off
%   octave> yyyy=2001
%   octave> [lon_value, lat_value, datestr_value, lai_value]=read_lai(yyyy);
%
% A file :file:`${PROJECT_ID}/LAI/laisen2001_write_lai_float.txt` can be written
% with data we just read::
%
%   octave> result=write_lai(lon_value, lat_value, datestr_value, lai_value);
%
% Files comparison might help detecting problems (either on files or in programs involved)::
%
%   $ diff ${PROJECT_ID}/LAI/laisen2001_float.txt ${PROJECT_ID}/LAI/laisen2001_write_lai_float.txt
%
% SEE ALSO
% ========
%
% :ref:`guide data LAI <data_lai>`
%
% :func:`simul_lai`
% :func:`write_lai`
%
% :func:`read_lai_2d`
%
% TODO
% ====
%
% status de retour
%
%
% description du résultat en ReST
%
% hard coded value of nb_data_line_max
% should be calculated according to latitude and longitude range and resolution
%
% check argument
%
% EVOLUTIONS
% ==========
%
% $Id: read_lai.m 325 2011-08-04 15:30:01Z pinsard $
%
% $URL$
%
% - fplod 20101202T141321Z aedon.locean-ipsl.upmc.fr (Darwin)
%
%   * new terminology and new format
%   * add parameter information in header
%
% - fplod 20100903T094128Z aedon.locean-ipsl.upmc.fr (Darwin)
%
%   * remove unique on lat and lon
%
% - fplod 20100827T103733Z aedon.locean-ipsl.upmc.fr (Darwin)
%
%   * creation to validate :file:`write_lai.m`
%
%-
%
global PROJECT_ID
%
% check argument+todo+
%
% build LAI ASCII filename to be read
fullfilename=[PROJECT_ID 'LAI/laisen'  num2str(yyyy) '_float.txt'];
disp(['iii : opening for reading ', fullfilename]);
% opening file
fid=fopen(fullfilename,'r');
%
% +todo+ check opening ok
%
% first line contains "Y,X" and an unknown quantity of days of years
% so we read this line and loop to find the number of columns
%
% get next line as a string
header_line = fgetl(fid);
%
% loop to find the number of columns
% and initialize the "day of year" array
index_column = 0;
% the maximum of doy is the number of days in year yyyy
% divided by 8 (according to time resolution of this dataset)
nb_doy_max=ceil((datenum(yyyy,12,eomday(yyyy,12)) - datenum(yyyy,1,1))/8);
doy=zeros(nb_doy_max,1);
doy(:)=NaN;
datestr_value=doy;
while ( ~isempty(header_line) )
   % parse next column label
   [next,header_line] = strtok(header_line,',');
   switch next
      case {'X'}
      case {'Y'}
      otherwise
         doy(index_column - 1) = str2double(next);
   end
   clear next
   % increase index of column
   index_column = index_column + 1;
end
%
% we know now how many columns are written on header
nb_column = index_column;
nb_doy = nb_column - 2;
clear index_column
%
% remove extra doy NaN values due to preallocation
doy=doy(~isnan(doy));
%
% initialisation d'un un tableau string de nb_doy elements à 9999999
datestr_value=repmat('99999999',nb_doy,1);
for index_doy=1:nb_doy
    datestr_value(index_doy,:)=datestr((datenum(yyyy,1,1) + doy(index_doy) - 1),'yyyymmdd');
end
%
% following lines contains lat, lon and LAI values
% we first preallocate arrays (boite,dx,dy formula or
% equivalent of IDL function file_lines ++)
% and initialize values to NaN
nb_data_line_max = 231002 - 2 ; %++
lat_value=zeros(nb_data_line_max,1);
lat_value(:)=NaN;
lon_value=zeros(nb_data_line_max,1);
lon_value(:)=NaN;
lai_line=zeros(nb_doy,nb_data_line_max);
lai_line(:,:)=NaN;
%
% we also need intermidiate arrays of lat and lon
lat_multi=zeros(nb_data_line_max,1);
lat_multi(:)=NaN;
lon_multi=zeros(nb_data_line_max,1);
lon_multi(:)=NaN;
%
% loop reading data from file (from second line till end of file)
index_data_line = 1;
while ( ~feof(fid) )
   % read one line at a time
   data_line = fgetl(fid)
   %
   % split data_line to values
   data=sscanf(data_line,'%f,');
   clear data_line
   % put values into arrays to be used for plot or computation
   lat_value(index_data_line)=data(1);
   lon_value(index_data_line)=data(2);
   lai_line(:,index_data_line)=data(3:nb_column);
   clear data
   %++disp(['data_line reading '  num2str(index_data_line)]);
   index_data_line = index_data_line + 1;
end
%
% close file
fclose(fid);
%
nb_data_line = index_data_line;
clear index_data_line
%
nb_lat=size(lat_value,1)
nb_lon=size(lon_value,1)
%
% organize LAI values in a 3D array (idate,ilat,ilon)
lai_value=zeros(nb_doy,nb_lat,nb_lon);
lai_value(:,:,:)=NaN;
for index_lat=1:nb_lat
   for index_lon=1:nb_lon
       index_data_line = index_lon + ((index_lat-1)* nb_lon);
       %++disp(['data_line storage '  num2str(index_data_line)]);
       lai_value(:,index_lat,index_lon)=lai_line(:,index_data_line);
   end
end
%
end