Matlab/TestMatlab

%======================================
% tests the flexpart matlab routines
%--------------------------------------
% HSO, 15.6.2007
%======================================

input='/nilu2/home/hso/data/flexpart/polarcat_test';

%reading data on/off
inp=1;

% setting display limits
lon=[-40 60];
lat=[0 89];
% convert to coordinates
lon1=lon+181;
lat1=lat+91;

% define a colormap
cmap=jet;
cmap(1,:)=[0.7 0.7 1.0];

% read the header
[h f]=flex_header(input);

% print out the header info
disp(h);

if f==1
  disp(['could not find header in ',input]);
else

  if inp==1
    % read one date
    data0=flex_read(h,0,0,15);
    % dimensions of data0:
    % data0(numxgrid, numygrid, numzgrid, nspec, nageclass, times);

    % calculate total column
    data1=flex_toc(h,data0);
    % dimensions of data1:
    % data1(numxgrid, numygrid, nspec, nageclass, times);

    % read several dates
    data2=flex_read(h,0,0,[3 7 21]);
    % calculate total column
    data3=flex_toc(h,data2);
    % dimensions of data3:
    % data3(numxgrid, numygrid, nspec, nageclass, times);

    % read all dates (slow!)
    %data4=flex_read(h);

  end % if inp

  %----------------------------------------------------------
  % plot the data
  %----------------------------------------------------------

  figure(1);
  set(gcf,'color','w');
  colormap(cmap);

  % first a snapshot from the single date input
  disp('Figure 1: a snapshot from the single date input');
  time=1;
  contourf(squeeze(data1(lon1(1):lon1(2),lat1(1):lat1(2),:,:,time))',50,'EdgeColor','none');
  set(gca,'XTick',1:10:100,'XTickLabel',lon(1):10:lon(2)');
  set(gca,'YTick',1:10:90,'YTickLabel',lat(1):10:lat(2)');
  title(sprintf('Total column concentration of species %s on %s',...
        h.species,datestr(h.dates(15))));
  xlabel('longitude');
  ylabel('latitude');

  figure(2)
  set(gcf,'color','w');
  colormap(cmap);

  % now two different dates from the multi-date input
  disp('Figure 2: two different dates from the multi-date input');
  subplot(1,2,1)
  time=1;
  contourf(squeeze(data3(lon1(1):lon1(2),lat1(1):lat1(2),:,:,time))',50,'EdgeColor','none');
  set(gca,'XTick',1:10:100,'XTickLabel',lon(1):10:lon(2)');
  set(gca,'YTick',1:10:90,'YTickLabel',lat(1):10:lat(2)');
  title(sprintf('Total column of species %s on %s',...
        h.species,datestr(h.dates(3))));
  xlabel('longitude');
  ylabel('latitude');

  subplot(1,2,2)
  time=3;
  contourf(squeeze(data3(lon1(1):lon1(2),lat1(1):lat1(2),:,:,time))',50,'EdgeColor','none');
  set(gca,'XTick',1:10:100,'XTickLabel',lon(1):10:lon(2)');
  set(gca,'YTick',1:10:90,'YTickLabel',lat(1):10:lat(2)');
  title(sprintf('Total column of species %s on %s',...
        h.species,datestr(h.dates(7))));
  xlabel('longitude');
  ylabel('latitude');

  disp('done');

end