/[MITgcm]/MITgcm_contrib/gael/matlab_class/gcmfaces_maps/m_map_gcmfaces.m
ViewVC logotype

Annotation of /MITgcm_contrib/gael/matlab_class/gcmfaces_maps/m_map_gcmfaces.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.8 - (hide annotations) (download)
Sat Jan 18 22:03:28 2014 UTC (11 years, 6 months ago) by gforget
Branch: MAIN
Changes since 1.7: +126 -28 lines
- add doCbar, doLabl, doFit optional parameters
- switch to Pacific centric maps for llc and cube
- add 1.1, 2.1, 3.1, 4.1 to 4.8 projection variants
- add hooks to xtrct
- add text diaply option
- add tag spec. to map & proj generated with this routine

1 gforget 1.1 function []=m_map_gcmfaces(fld,varargin);
2 gforget 1.3 %object: gcmfaces front end to m_map
3 gforget 1.7 %inputs: fld is the 2D field to be mapped, or a cell (see below).
4 gforget 1.4 %optional: proj is either the index (integer; 0 by default) of pre-defined
5     % projection(s) or parameters to pass to m_proj (cell)
6 gforget 1.7 %more so: other optional paramaters can be provided ONLY AFTER proj,
7 gforget 1.4 % and must take the following form {'name',param1,param2,...}
8 gforget 1.6 % those that are currently active are
9 gforget 1.4 % {'myCaxis',myCaxis} is the color axis ('auto' by default)
10     % {'do_m_coast',do_m_coast} adds call to m_coast (1; default) or not (0).
11 gforget 1.6 % {'myCmap',myCmap} is the colormap name ('jet' by default)
12 gforget 1.7 % {'doHold',1} indicates to 'hold on' and superimpose e.g. a contour
13 gforget 1.8 % {'doCbar',1} indicates to include the colorbar (1 by default)
14     % {'doLabl',1} indicates to include the axes labels (1 by default)
15     % {'doFit',1} indicates to exclude white space (0 by default)
16 gforget 1.3 %
17 gforget 1.4 %notes: - for proj==0 (i.e. the default) three panels will be plotted :
18     % a mercator projection over mid-latitude, and two polar
19 gforget 1.7 % stereographic projections. The other predefined proj are
20 gforget 1.3 % -1 mollweide or cylindrical projection
21     % 1 mercator projection only
22     % 2 arctic stereographic projection only
23     % 3 antarctic stereographic projection only
24 gforget 1.7 % 1.1 and 1.2 are atlantic mercator projections
25 gforget 1.4 % - myTitle is currently not used; it will take the form
26     % {'myTitle',myTitle} is the title (none by default)
27 gforget 1.7 % - myCaxis can be specified with more than 2 values, in which
28 gforget 1.8 % case gcmfaces_cb will be used.
29 gforget 1.7 % - if fld is a 2D field we use pcolor to display it
30     % - if fld is a cell one can specify the plotting tool.
31     % The cell must then have the form {plotType,FLD,varargin}
32 gforget 1.8 % where plotType is e.g. 'pcolor' or 'contour', FLD is
33     % the 2D field to be plotted, and varargin are options
34 gforget 1.7 % to pass over to the contour command.
35 gforget 1.8 % - Hence .e.g. fld={'contour',mygrid.Depth,'r'} will draw
36 gforget 1.7 % red contours, while fld=mygrid.Depth will color shade.
37 gforget 1.1
38 gforget 1.2 %check that m_map is in the path
39     aa=which('m_proj'); if isempty(aa); error('this function requires m_map that is missing'); end;
40    
41 gforget 1.1 global mygrid;
42    
43 gforget 1.4 %get optional parameters
44     if nargin>1; proj=varargin{1}; else; proj=0; end;
45 gforget 1.7 if iscell(proj);
46     error('not yet implemented');
47 gforget 1.4 else;
48     choicePlot=proj;
49     end;
50 gforget 1.7 %determine the type of plot
51     if iscell(fld); myPlot=fld{1}; else; myPlot='pcolor'; fld={'pcolor',fld}; end;
52 gforget 1.4 %set more optional paramaters to default values
53 gforget 1.8 myCaxis=[]; myTitle=''; do_m_coast=1; myCmap='jet'; doHold=0; doCbar=1; doLabl=1; doFit=0;
54 gforget 1.4 %set more optional paramaters to user defined values
55     for ii=2:nargin-1;
56     if ~iscell(varargin{ii});
57     warning('inputCheck:m_map_gcmfaces_1',...
58     ['As of june 2011, m_map_gcmfaces expects \n'...
59     ' its optional parameters as cell arrays. \n'...
60     ' Argument no. ' num2str(ii+1) ' was ignored \n'...
61     ' Type ''help m_map_gcmfaces'' for details.']);
62     elseif ~ischar(varargin{ii}{1});
63     warning('inputCheck:m_map_gcmfaces_2',...
64     ['As of june 2011, m_map_gcmfaces expects \n'...
65     ' its optional parameters cell arrays \n'...
66     ' to start with character string. \n'...
67     ' Argument no. ' num2str(ii+1) ' was ignored \n'...
68     ' Type ''help m_map_gcmfaces'' for details.']);
69     else;
70     if strcmp(varargin{ii}{1},'myCaxis')|...
71 gforget 1.6 strcmp(varargin{ii}{1},'myCmap')|...
72 gforget 1.7 strcmp(varargin{ii}{1},'doHold')|...
73 gforget 1.8 strcmp(varargin{ii}{1},'doCbar')|...
74     strcmp(varargin{ii}{1},'doLabl')|...
75     strcmp(varargin{ii}{1},'doFit')|...
76 gforget 1.4 strcmp(varargin{ii}{1},'do_m_coast')|...
77     strcmp(varargin{ii}{1},'myTitle');
78     eval([varargin{ii}{1} '=varargin{ii}{2};']);
79     else;
80     warning('inputCheck:m_map_gcmfaces_3',...
81     ['unknown option ''' varargin{ii}{1} ''' was ignored']);
82     end;
83     end;
84     end;
85 gforget 1.7
86 gforget 1.4 %make parameter inferences
87 gforget 1.7 if length(myCaxis)==0;
88     plotCBAR=0;
89     elseif length(myCaxis)==2;
90     plotCBAR=1;
91     else;
92     plotCBAR=2;
93     end;
94     %
95     if choicePlot==0&~doHold;
96     clf;
97     elseif ~doHold;
98     cla;
99     else;
100     hold on;
101     end;
102    
103     %re-group param:
104     param.plotCBAR=plotCBAR;
105     param.myCmap=myCmap;
106     param.myCaxis=myCaxis;
107     param.do_m_coast=do_m_coast;
108     param.doHold=doHold;
109 gforget 1.8 param.doCbar=doCbar;
110     param.doLabl=doLabl;
111     param.doFit=doFit;
112 gforget 1.7 param.myPlot=myPlot;
113    
114     %do the plotting:
115     if (choicePlot~=0&choicePlot~=1&choicePlot~=2&choicePlot~=3);
116     do_my_plot(fld,param,choicePlot);
117     end;%if choicePlot==0|choicePlot==1;
118 gforget 1.1
119     if choicePlot==0; subplot(2,1,1); end;
120 gforget 1.7 if choicePlot==0|choicePlot==1;
121 gforget 1.8 if mygrid.nFaces~=5&mygrid.nFaces~=6;
122     do_my_plot(fld,param,1);
123     else;
124     do_my_plot(fld,param,1.1);
125     end;
126 gforget 1.7 end;%if choicePlot==0|choicePlot==1;
127 gforget 1.3
128 gforget 1.1 if choicePlot==0; subplot(2,2,3); end;
129 gforget 1.7 if choicePlot==0|choicePlot==2;
130     do_my_plot(fld,param,2);
131     end;%if choicePlot==0|choicePlot==1;
132 gforget 1.1
133     if choicePlot==0; subplot(2,2,4); end;
134 gforget 1.7 if choicePlot==0|choicePlot==3;
135     do_my_plot(fld,param,3);
136     end;%if choicePlot==0|choicePlot==1;
137 gforget 1.1
138 gforget 1.8 if plotCBAR==2&strcmp(myPlot,'pcolor')&doCbar;
139 gforget 1.6 cbar=gcmfaces_cmap_cbar(myCaxis,{'myCmap',myCmap});
140 gforget 1.7 if choicePlot==0;
141 gforget 1.3 set(cbar,'Position',[0.92 0.15 0.02 0.75]);
142     elseif choicePlot==-1;
143     set(cbar,'Position',[0.92 0.35 0.02 0.3]);
144 gforget 1.7 elseif choicePlot==1;
145 gforget 1.3 set(cbar,'Position',[0.92 0.34 0.02 0.35]);
146     else;
147     set(cbar,'Position',[0.92 0.23 0.02 0.55]);
148 gforget 1.7 end;
149     end;
150    
151 gforget 1.8 if doFit;
152     if doLabl&~doCbar;
153     set(gca,'LooseInset',[0.05 0.02 0.03 0.03]);
154     else;
155     set(gca,'LooseInset',[0.01 0 0.03 0.03]);
156     end;
157     tmp1=get(gca,'PlotBoxAspectRatio');
158     set(gcf,'PaperPosition',[0 4 8 8/tmp1(1)]);
159     end;
160 gforget 1.7
161     function []=do_my_plot(fld,param,proj);
162    
163     gcmfaces_global;
164    
165 gforget 1.8 %default m_grid params:
166     xloc='bottom'; xtic=[-180:60:180]; xticlab=1;
167     yloc='left'; ytic=[-60:30:60]; yticlab=1;
168    
169     %choice of projection:
170 gforget 1.7 if proj==-1;
171     %%m_proj('Miller Cylindrical','lat',[-90 90]);
172     %m_proj('Equidistant cylindrical','lat',[-90 90]);
173     %m_proj('mollweide','lon',[-180 180],'lat',[-80 80]);
174     m_proj('mollweide','lon',[-180 180],'lat',[-88 88]);
175 gforget 1.8 myConv='pcol'; xticlab=0; yticlab=0;
176 gforget 1.7 elseif proj==1;
177     m_proj('Mercator','lat',[-70 70]);
178 gforget 1.8 myConv='pcol';
179     elseif proj==1.1;
180     m_proj('Mercator','lat',[-70 70],'lon',[0 360]+20);
181     myConv='pcol';
182     xtic=[-360:60:360]; ytic=[-60:30:60];
183 gforget 1.7 elseif proj==2;
184     m_proj('Stereographic','lon',0,'lat',90,'rad',40);
185 gforget 1.8 myConv='convert2arctic';
186     yloc='bottom'; ytic=[50:10:90];
187     elseif proj==2.1;
188     m_proj('Stereographic','lon',0,'lat',90,'rad',60);
189     myConv='convert2arctic';
190     yloc='bottom'; ytic=[30:10:90];
191 gforget 1.7 elseif proj==3;
192     m_proj('Stereographic','lon',0,'lat',-90,'rad',40);
193 gforget 1.8 myConv='convert2southern';
194     xloc='top'; ytic=[-90:10:-50];
195     elseif proj==3.1;
196     m_proj('Stereographic','lon',0,'lat',-90,'rad',60);
197     myConv='convert2southern';
198     xloc='top'; ytic=[-90:10:-30];
199     elseif proj==4.1;
200     m_proj('mollweide','lat',[25 75],'lon',[-100 30]);
201     myConv='pcol';
202     xtic=[-100:20:30]; ytic=[30:10:70];
203     elseif proj==4.2;
204     m_proj('mollweide','lat',[-30 30],'lon',[-65 20]);
205     myConv='pcol';
206     xtic=[-60:10:20]; ytic=[-30:10:30];
207     elseif proj==4.3;
208     m_proj('mollweide','lat',[-75 -25],'lon',[-75 25]);
209     myConv='pcol';
210     xtic=[-70:20:20]; ytic=[-70:10:-30]; xloc='top';
211     elseif proj==4.4;
212     m_proj('mollweide','lat',[25 75],'lon',[-240 -110]);
213     myConv='pcol';
214     xtic=[-240:20:-120]; ytic=[30:10:60];
215     elseif proj==4.5;
216     m_proj('mollweide','lat',[-30 30],'lon',[-240 -70]);
217     myConv='pcol';
218     xtic=[-240:20:-70]; ytic=[-30:10:30];
219     elseif proj==4.6;
220     m_proj('mollweide','lat',[-75 -25],'lon',[-215 -60]);
221     myConv='pcol';
222     xtic=[-210:20:-70]; ytic=[-70:10:-30]; xloc='top';
223     elseif proj==4.7;
224     m_proj('mollweide','lat',[-30 30],'lon',[15 155]);
225     myConv='pcol';
226     xtic=[10:20:160]; ytic=[-30:10:30]; xloc='top';
227     elseif proj==4.8;
228     m_proj('mollweide','lat',[-75 -25],'lon',[15 155]);
229     myConv='pcol';
230     xtic=[10:20:160]; ytic=[-70:10:-30]; xloc='top';
231     else;
232     error('undefined projection');
233 gforget 1.3 end;
234 gforget 1.1
235 gforget 1.8 if ~param.doLabl; xticlab=0; yticlab=0;end; %omit labels
236    
237     m_grid_opt=['''XaxisLocation'',xloc,''YaxisLocation'',yloc,''xtick'',xtic,''ytick'',ytic'];
238     if xticlab==0; m_grid_opt=[m_grid_opt ',''xticklabel'',[]']; end;
239     if yticlab==0; m_grid_opt=[m_grid_opt ',''yticklabel'',[]']; end;
240    
241 gforget 1.7 if strcmp(param.myPlot,'pcolor')|strcmp(param.myPlot,'contour');
242 gforget 1.8 x=mygrid.XC;
243     %mask out the XC padded zeros
244     if isfield(mygrid,'xtrct');
245     pt1=mygrid.xtrct.pt1face; pt2=mygrid.xtrct.pt2face;
246     if pt1~=pt2;
247     for iF=1:6; if iF~=pt1&iF~=pt2; x{iF}(:)=NaN; end; end;
248     end;
249     x(x<0)=x(x<0)+360;
250     end;
251     if strcmp(myConv,'pcol');
252     [xx,yy,z]=convert2pcol(x,mygrid.YC,fld{2});
253 gforget 1.7 else;
254 gforget 1.8 eval(['xx=' myConv '(x);']);
255     eval(['yy=' myConv '(mygrid.YC);']);
256     eval(['z=' myConv '(fld{2});']);
257 gforget 1.7 end;
258     [x,y]=m_ll2xy(xx,yy);
259     if strcmp(param.myPlot,'pcolor');
260     if sum(~isnan(x(:)))>0; pcolor(x,y,z); shading interp; end;
261 gforget 1.8 if param.plotCBAR==0;
262     colormap(param.myCmap); if param.doCbar; colorbar; end;
263     elseif param.plotCBAR==1;
264     caxis(param.myCaxis); colormap(param.myCmap); if param.doCbar; colorbar; end;
265     else;
266     cbar=gcmfaces_cmap_cbar(param.myCaxis,{'myCmap',param.myCmap}); delete(cbar);
267 gforget 1.7 end;
268     if param.do_m_coast; m_coast('patch',[1 1 1]*.7,'edgecolor','none'); end;
269 gforget 1.8 eval(['m_grid(' m_grid_opt ');']);
270 gforget 1.7 elseif strcmp(param.myPlot,'contour');
271     if ~param.doHold;
272     if param.do_m_coast; m_coast('patch',[1 1 1]*.7,'edgecolor','none'); end;
273 gforget 1.8 eval(['m_grid(' m_grid_opt ');']);
274 gforget 1.7 end;
275     if length(fld)==2; fld{3}='k'; end;
276     hold on; contour(x,y,z,fld{3:end});
277     end;
278     elseif strcmp(param.myPlot,'plot');
279     if ~param.doHold;
280     if param.do_m_coast; m_coast('patch',[1 1 1]*.7,'edgecolor','none'); end;
281 gforget 1.8 eval(['m_grid(' m_grid_opt ');']);
282 gforget 1.7 end;
283     [x,y]=m_ll2xy(fld{2},fld{3});
284     hold on; plot(x,y,fld{4:end});
285 gforget 1.8 elseif strcmp(param.myPlot,'text');
286     if ~param.doHold;
287     if param.do_m_coast; m_coast('patch',[1 1 1]*.7,'edgecolor','none'); end;
288     eval(['m_grid(' m_grid_opt ');']);
289     end;
290     [x,y]=m_ll2xy(fld{2},fld{3});
291     if length(fld)>4; cc=fld{5}; else; cc='k'; end;
292     hold on; hold on; text(x,y,fld{4},'Color',cc,'FontSize',16,'FontWeight','bold');
293 gforget 1.7 end;
294 gforget 1.1
295 gforget 1.8 %add tag spec. to map & proj generated with this routine
296     set(gca,'Tag',['gfMap' num2str(proj)]);
297    

  ViewVC Help
Powered by ViewVC 1.1.22