[Dart-dev] [3845] DART/trunk/models/cam/matlab/DART_xcorrel.m: This script performs the correlation of one variable at a particular

nancy at ucar.edu nancy at ucar.edu
Mon May 4 17:00:13 MDT 2009


An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/dart-dev/attachments/20090504/cf614e83/attachment-0001.html 
-------------- next part --------------
Added: DART/trunk/models/cam/matlab/DART_xcorrel.m
===================================================================
--- DART/trunk/models/cam/matlab/DART_xcorrel.m	                        (rev 0)
+++ DART/trunk/models/cam/matlab/DART_xcorrel.m	2009-05-04 23:00:13 UTC (rev 3845)
@@ -0,0 +1,160 @@
+function DART_xcorrel(datadir,fname,varname1,location,varname2,level2)
+% DART_xcorrel  explores correlation between one DART location and the DART ensemble
+%
+% datadir  = '/ptmp/raeder/Cam3.6/Hist0';
+% fname    = 'Pr_06_ens.nc';
+% varname1 = 'T';
+% location = [284, -54, 600];  % lon/lat/level
+% varname2 = 'T';
+% level2   = 600;
+%
+% DART_xcorrel(datadir,fname,varname1,location,varname2,level2)
+
+
+% Data Assimilation Research Testbed -- DART
+% Copyright 2004-2007, Data Assimilation Research Section
+% University Corporation for Atmospheric Research
+% Licensed under the GPL -- www.gpl.org/licenses/gpl.html
+%
+% <next few lines under version control, do not edit>
+% $URL$
+% $Id$
+% $Revision$
+% $Date$
+
+
+posmat = [0.1 0.60 0.8 0.25;
+          0.1 0.10 0.8 0.40];
+
+%----------------------------------------------------------------------
+% Get the "time series" from the DART file:
+%----------------------------------------------------------------------
+
+% For this variable, should query for coordinate variables
+% not handling variables on staggered grids correctly as it stands.
+
+camname = sprintf('%s/%s',datadir,fname);
+
+CAMlevels = nc_varget(camname,'lev'); CAMNlevels = length(CAMlevels);
+CAMlats   = nc_varget(camname,'lat'); CAMNlats   = length(CAMlats);
+CAMlons   = nc_varget(camname,'lon'); CAMNlons   = length(CAMlons);
+CAMunits  = nc_attget(camname,varname1,'units'); 
+CAMtimestr = timebase(camname);
+
+diffs  = abs(CAMlevels - location(3)); [camlevel,camlevind] = min(diffs);
+diffs  = abs(CAMlats   - location(2)); [camlat,  camlatind] = min(diffs);
+diffs  = abs(CAMlons   - location(1)); [camlon,  camlonind] = min(diffs);
+
+camlevel = CAMlevels(camlevind);
+camlat   = CAMlats(  camlatind);
+camlon   = CAMlons(  camlonind);
+
+%   float T(time, copy, lat, lon, lev) ;
+%                T:long_name = "Temperature" ;
+%                T:units = "K" ;
+
+start = [ 0   0   0   0   camlevind-1];
+count = [ 1   -1  -1  -1        1];
+
+bob = nc_varget(camname,varname1,start,count);
+x   = bob(:,camlatind,camlonind); 
+
+cammean = squeeze(mean(bob,1));
+
+%----------------------------------------------------------------------
+% Get the DART fields.
+% float Q(time, copy, lat, lon, lev) ;
+%----------------------------------------------------------------------
+
+prfname    = sprintf('%s/%s',datadir,fname);
+DARTlons   = nc_varget(prfname,'lon');
+DARTlats   = nc_varget(prfname,'lat');
+DARTlevels = nc_varget(prfname,'lev');
+
+diffs = abs(DARTlevels - level2);
+[dartlevel,dartlevind] = min(diffs);
+dartlevel = DARTlevels(dartlevind);
+
+start = [ 0  0  0  0 dartlevind-1];
+count = [-1 -1 -1 -1      1];
+y     = nc_varget(prfname,varname2,start,count);
+
+[Nmem, Nlat, Nlon] = size(y);
+
+%----------------------------------------------------------------------
+% Reshape the DART data to be in state-space form ... Nt-by-Nx
+% perform the correlation,
+% reshape the correlation coefficients back to a Lat-Lon matrix.
+%----------------------------------------------------------------------
+
+Yss = reshape(y,[Nmem Nlat*Nlon]); 
+a   = corr(x,Yss);
+b   = reshape(a,[Nlat Nlon]);
+
+clear Yss a
+
+bob = b(:);
+abob = abs(bob);
+datamax = max(abob);
+
+%datamax = 0.5;
+%----------------------------------------------------------------------
+% make the plots
+%----------------------------------------------------------------------
+
+figure(1); clf; orient tall
+
+str1 = sprintf('%s correlated against %s',varname1,varname2);
+str2 = sprintf('%s [%.2f , %.2f] level %.2f',varname1,camlon,camlat,camlevel);
+str3 = sprintf('%s level %.2f',varname2,dartlevel);
+
+subplot('position',posmat(1,:))
+plot(x)
+title( {str2, CAMtimestr} )
+ylabel(CAMunits)
+xlabel('Ensemble Member');
+
+subplot('position',posmat(2,:))
+mymap = jet;
+colormap(mymap)
+imagesc(CAMlons,CAMlats,cammean)
+set(gca,'YDir','normal')
+h = colorbar;
+set(get(h,'YLabel'),'String',CAMunits)
+worldmap;
+axis image
+hold on; 
+plot([  0 360],[camlat camlat],'k', [camlon camlon],[-90 90],'k');
+title(sprintf('Ensemble mean %s level %.2f   %s',varname1,camlevel,CAMtimestr))
+
+
+figure(2); clf; orient tall
+
+
+subplot('position',posmat(2,:))
+imagesc(DARTlons,DARTlats,b,[-datamax datamax])
+set(gca,'YDir','normal')
+mymap = jet;
+mymap(26:39,:) = 1;
+colormap(mymap);
+colorbar
+worldmap;
+axis image
+hold on; 
+plot([  0 360],[camlat camlat],'k', [camlon camlon],[-90 90],'k');
+
+title({str1,str2,str3,CAMtimestr})
+
+
+
+
+function timestring = timebase(fname)
+
+t      = nc_varget(fname,'time');
+units  = nc_attget(fname,'time','units'); 
+mycal  = nc_attget(fname,'time','calendar'); 
+
+timebase = sscanf(units,'%*s%*s%d%*c%d%*c%d'); % YYYY MM DD
+timeorigin = datenum(timebase(1),timebase(2),timebase(3));
+timestring = datestr(t + timeorigin);
+


Property changes on: DART/trunk/models/cam/matlab/DART_xcorrel.m
___________________________________________________________________
Name: svn:mime-type
   + text/x-matlab
Name: svn:keywords
   + Date Revision Author HeadURL Id
Name: svn:eol-style
   + native


More information about the Dart-dev mailing list