source: trunk/src/scatterplot.m @ 325

function [longi,latit,amsutot]=scatterplot(nosat, mois,jour,heure,lonmin,lonmax,latmin,latmax)

%SCATTERPLOT +todo+
%+
%
% =============
% scatterplot.m
% =============
%
% .. function:: scatterplot(mois,jour,heure,lonmin,lonmax,latmin,latmax)
%
% DESCRIPTION
% ===========
%
%   :param mois: mois 
%   :type mois: array of int8
%   :raise mois: required
%
%   :param jour: jour debut
%   :type jour: array of int8
%   :raise jour: required
%
%   :param heure: heure ++ amsu
%   :type heure: int8
%   :raise heure: required
%
%   :param lonmin: longitude min, W < 0
%   :units lonmin: deg
%   :type lonmin: double
%   :raise lonmin: required
%
%   :param lonmax: longitude max, W < 0
%   :units lonmax: deg
%   :type lonmax: double
%   :raise lonmax: required
%
%   :param latmin: latitude min, N > 0
%   :units latmin: deg
%   :type latmin: double
%   :raise latmin: required
%
%   :param latmax: latitude max, N > 0
%   :units latmax: deg
%   :type latmax: double
%   :raise latmax: required
%
% tentative de cumul de quelques fichiers AMSU et MSG pour traces
% histogrammes et diagramme de dispersion
%
% +todo+ lit
%
% +todo+ calcule
%
% +todo+ dessine
%
% EXAMPLES
% ========
%
% To plot +todo+::
%
%  octave> clear all
%  octave> varamma_startup
%  octave> mois=[8 8 8 8 8 8];
%  octave> jour=[7 18 31 25 22 02];
%  octave> heure=[191517 181344 80216 184655 181839 175432];
%  octave> lonmin=-20;
%  octave> lonmax=-10;
%  octave> latmin=10;
%  octave> latmax=20;
%  octave> [longi,latit,amsutot]=scatterplot(mois,jour,heure,lonmin,lonmax,latmin,latmax);
%
% SEE ALSO
% ========
%
% :ref:`guide data MSG <data_msg>`
% :ref:`guide data AMSU <data_amsu>`
%
% :ref:`varamma_startup.m`
%
% TODO
% ====
%
% get rid of mlint output
%
% yyyy arg
%
%
% status de retour
%
% split reading and calculate and plot
%
% check argument
%
% EVOLUTIONS
% ==========
%
% $Id$
%
% $URL$
%
% - fplod 20110121T140728Z aedon.locean-ipsl.upmc.fr (Darwin)
%
%   * complement on header
%
% - fplod 20100901T150719Z aedon.locean-ipsl.upmc.fr (Darwin)
%   
%   * restore distinction between ASMU and MSG directories
%   * replace int2str by num2str to be able to specify a format
%   * revision of indentation
%
% - fplod 20100604T092417Z aedon.locean-ipsl.upmc.fr (Darwin)
%
%   * correction for running with matlab
%   * minimize mlint output
%
% - fplod 20100604T083004Z aedon.locean-ipsl.upmc.fr (Darwin)
%
%   * add usage of varamma_startup : no more chemin parameter, PROJECT_ID value
%     used instead
%   * no more distinction between ASMU and MSG directories : each file is
%     intended to be in PROJECT_ID (links are welcome !)
%
% - fplod 20100521T150917Z aedon.locean-ipsl.upmc.fr (Darwin)
%
%   * minimal header
%
% - lelod 20100521
%
%   * creation
%-

global PROJECT_ID;
global PROJECT_OD;

nbfic=size(jour(:),1);
numheure=fix(heure(:)/10000.);
minute=fix((heure(:)-numheure(:)*10000)/100)/4;
quart=minute(:)/4;
ih=fix(4*numheure(:)+1+quart(:));

% AMSU
z_deb=[latmin lonmin];
z_fin=[latmax lonmax];
amsutot=[];
yyyy=2006;
for i=1:nbfic
   mo=num2str(mois(i),'%2.2d');
   jo=num2str(jour(i),'%2.2d');
   he=num2str(heure(i),'%6.6d');
   sat=num2str(nosat(i),'%2.2d');
   file=['extracted-amsubch5-noaa',sat,'_', num2str(yyyy,'%4.4d'), mo, jo,'-',he,'_afratl_Nadir.nc'];
%   file=[PROJECT_ID,'AMSU/2006/tmp/', num2str(yyyy,'%4.4d'), '/',
%   mo, '/', file];
  file=[PROJECT_ID,'AMSU/tmp/', file];

   file=netcdf(file,'nowrite');
   time=file{'time'};
   lon=file{'lon'};
   lat=file{'lat'};
   tbamsu=file{'tbch5b'};
   mask=file{'mask'};
   la=find(lat(:)>z_deb(1) & lat(:)<z_fin(1));
   lo=find(lon(:)>z_deb(2) & lon(:)<z_fin(2));
   latit=lat(la);
   longi=lon(lo);
   tbzone=squeeze(tbamsu(1,la,lo));
   tbzone(find(tbzone == -999)) = NaN;
   if i==1
     amsutot=tbzone;
   end
   if i>=2
     amsutot=[amsutot;tbzone];
   end
end

% MSG
msgtot=[];
for i=1:nbfic
   mo=num2str(mois(i),'%2.2d');
   jo=num2str(jour(i),'%2.2d');
   date=[num2str(yyyy,'%4.4d'),'-',mo,'-',jo];
   file=['extracted-msg-tb108_',date,'_15min.nc'];
   file=[PROJECT_ID,'MSG/', num2str(yyyy,'%4.4d'), '/', mo, '/', file];
   file=netcdf(file,'nowrite');
   tbmsg=file{'DATA'};
   lon=file{'LONGITUDE'};
   lat=file{'LATITUDE'};
   time=file{'TIME'};
   la=find(lat(:)>z_deb(1) & lat(:)<z_fin(1));
   lo=find(lon(:)>z_deb(2) & lon(:)<z_fin(2));
   latit=lat(la);
   longi=lon(lo);
   hh=ih(i);
   tbzon=0.01*squeeze(tbmsg(hh,la,lo))+327.68;
   if i==1
     msgtot=tbzon;
   end
   if i>=2
     msgtot=[msgtot;tbzon];
   end
end


%trace histogrammes
figure(2)
Tseuil1=280;   % on enleve les points chauds (surface)
a=find(msgtot <Tseuil1);
hist(msgtot(a))
xlabel('TbMSG');
figure(3)
Tseuil2=250;   % on enleve les points chauds (surface)
a=find(amsutot <Tseuil2);
hist(amsutot(a))
xlabel('TbAMSU');

%trace diagramme de dispersion tbamsu vs tbmsg
figure(1)
a=find((msgtot <Tseuil1) & (amsutot <Tseuil2));
%scatter(amsutot(a),msgtot(a), "+")
plot(amsutot(a),msgtot(a), '+', 'markersize', 1);
xlabel('TbAMSU');
ylabel('TbMSG');