function [X,Y,NUMPOINTS,RADIUS,MEANTEMP,STDTEMP,MEANSAL,STDSAL,MEANST,STDST,MEANC,STDC]=cl_climo_planview(filename,season,upperdepth,lowerdepth,gridspacing,Ro,minnumpts,radiusinc,cutoffhighdepth,cutofflowdepth,westbound,eastbound,southbound,northbound,qc,plotyes); % % CL_CLIMO_PLANVIEW This program takes a file of unevenly spaced hydrographic % data and creates a planview climatology based on user-specified inputs. % % function [X,Y,NUMPOINTS,RADIUS,MEANTEMP,STDTEMP,MEANSAL,STDSAL,MEANST,STDST,MEANC,STDC]=cl_climo_planview(filename,season,upperdepth,lowerdepth,gridspacing,Ro,minnumpts,radiusinc,cutoffhighdepth,cutofflowdepth,westbound,eastbound,southbound,northbound,qc,plotyes); % % Supporting m-files: % - Matlab mapping toolbox % - The m_map toolbox, free at this website: http://www2.ocgy.ubc.ca/~rich/map.html % - The CSIRO Matlab Seawater library, free at this website: % http://www.marine.csiro.au/datacentre/ext_docs/CSIRO_MatLAB_Seawater_Library.pdf % % Inputs: % filename - specifies name of data file % Data file, in ASCII text format, must have the following columns: % 1. Decimal longitude % 2. Decimal latitude % 3. Year % 4. Month % 5. Depth (meters) % 6. Temperature (C) % 7. Salinity (PSU) % season string - is a listing of all seasons (by number) to include in average. % For example, Jan & Feb would be '1 2'. Enter 'All' for all. % upperdepth upper bound (meters) of the climatology. % lowerdepth maximum depth (meters) of the climatology. % gridspacing grid spacing in degrees (1/6th degree is 19km lat / 14km lon resolution at 40N). % Ro Search radius in km (typical is 18km for 0.2 degree grid spacing). % In this version of the program the search radius will % increase iteratively (increment of 10km) until the minimum % minnumpts is met. % minnumpts minimum number of data points to use for an average. If % this number is not met, search radius will increase. % radiusinc Increment that radius will be increased in km to meet minnumpts % If 0, radius will not be iteratively increased. Output nodes % that don't have minnumpts observations will be given NaNs. % cutoffhighdepth Cutoff depth (shallow) - stations with a depth SHALLOWER than this value will be excluded. % cutofflowdepth Cutoff depth (deep) - stations with a depth DEEPER than this value will be excluded. % westbound western domain bound - decimal degrees. % eastbound eastern domain bound - decimal degrees. % southbound southern domain bound - decimal degrees. % northbound northern domain bound - decimal degrees. % qc 0 or 1; 0 will skip the quality control routine and 1 will run it. % plotyes 0, 1, or 2; 0 will suppress plots, 1 will show plots, and 2 will show plots and save them out as PNG files to the current directory. % % Outputs: % X Cross-shelf distance of output grid points in KM (meshgrid format). % Y Depth of output grid points in M (meshgrid format). % NUMPOINTS Number of points per output grid node. This, like all of the other output fields, are in the meshgrid format for easy plotting. For example, pcolor(X,Y,NUMPOINTS). % RADIUS Final search radius in KM. % MEANTEMP Mean temperature (C). % STDTEMP Standard deviation of temperature. % MEANSAL Mean salinity (PSU). % STDSAL Standard deviation of salinity. % MEANST Mean sigma-theta (density) (kg/m3) (sigma-theta and sound speed outputs computed using a 1024-realization monte-carlo routine with mean & std temp/salinity as inputs). % STDST Standard deviation of sigma-theta. % MEANC Mean sound speed (m/s). % STDC Standard deviation of sound speed. % % Example: Nantucket Shoals test case % [X,Y,NUMPOINTS,RADIUS,MEANTEMP,STDTEMP,MEANSAL,STDSAL,MEANST,STDST,MEANC,STDC]=cl_climo_planview('ns-all-data.dat','8 9',0,15,.2,18,10,10,50,2500,-72,-69,39,41,1,1); % % Given a data file 'ns-all-data.dat', a climatology will be created for August % and September near-surface (0 to 15 meters deep) data. The output grid % will have 0.2 degree spacing and a search radius of 18km. If the output % node doesn't have at least 10 realizations, the search radius will be % increased in 10km increments until 10 data points are included. Stations % shallower than 50m or deeper than 2500m will be excluded. The study area % bounds are 72-69 West, 39-41 North. The quality control routine will be % run and the data will be plotted but not saved. % % Version 5.0, 8Nov2006, Chris Linder clinder@whoi.edu % http://www.whoi.edu/science/PO/people/clinder/software_climo_planview.html % % This program is free software; you can redistribute it and/or % modify it under the terms of the GNU General Public License % as published by the Free Software Foundation; either version 2 % of the License, or (at your option) any later version. % % This program is distributed in the hope that it will be useful, % but WITHOUT ANY WARRANTY; without even the implied warranty of % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the % GNU General Public License for more details. more off; if plotyes==2 saveprefix=[filename(1:end-4),'_']; end; disp(' '); eval(sprintf('disp(''Loading %s...'');',filename)); pause(0.1); eval(sprintf('rawdata=load(''%s'');',filename)); disp(' '); disp('Extracting casts... '); % Extract 5 min Terrainbase bathymetry file for this area (m_map toolbox) m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); [data,LON,LAT]=m_tbase([westbound eastbound southbound northbound]); if season(1)=='A' numseasons=[1:1:12]; else numseasons=str2num(season); end; i1=find(rawdata(:,4)==1); i2=find(rawdata(:,4)==2); i3=find(rawdata(:,4)==3); i4=find(rawdata(:,4)==4); i5=find(rawdata(:,4)==5); i6=find(rawdata(:,4)==6); i7=find(rawdata(:,4)==7); i8=find(rawdata(:,4)==8); i9=find(rawdata(:,4)==9); i10=find(rawdata(:,4)==10); i11=find(rawdata(:,4)==11); i12=find(rawdata(:,4)==12); ichoice=[]; for sealoop=1:size(numseasons,2) eval(sprintf('ichoice=[ichoice; i%d];',numseasons(sealoop))); end; % Setup output grid dimensions newlat=northbound : -1 * abs(gridspacing) : southbound; newlon=westbound : abs(gridspacing) : eastbound; [X,Y]=meshgrid(newlon,newlat); seadist=[]; ALLPOINTS=rawdata(ichoice,:); placeholder=1; allcasts=[]; for rownum=2:size(ALLPOINTS,1) % 1 possibility - lat OR lon OR month OR year changes OR depth decreases: % then it's a new cast... if ALLPOINTS(rownum,1) ~= ALLPOINTS(rownum-1,1) | ALLPOINTS(rownum,2) ~= ALLPOINTS(rownum-1,2) | ALLPOINTS(rownum,3) ~= ALLPOINTS(rownum-1,3) | ALLPOINTS(rownum,4) ~= ALLPOINTS(rownum-1,4) | ALLPOINTS(rownum,5) < ALLPOINTS(rownum-1,5) castprops=ALLPOINTS(placeholder:rownum-1,:); % Need to check for duplicates in the depth record (rare) and % eliminate them [b,i,j]=unique(castprops(:,5)); castprops=castprops(i,:); % Compute profiles of sound speed and sigma-t and add to castprops % array castprops=[castprops sw_svel(castprops(:,7),castprops(:,6),abs(castprops(:,5))) sw_dens(castprops(:,7),castprops(:,6),abs(castprops(:,5)))]; % Find cast water depth depthofprofile=interp2(LON,LAT,data,castprops(1,1),castprops(1,2)); if depthofprofile >= -abs(cutofflowdepth) & depthofprofile <= -abs(cutoffhighdepth) & size(castprops,1) > 1 % Cast must be within shallow/deep bounds and have more than one point! i_inrange=find(abs(castprops(:,5)) <= lowerdepth & abs(castprops(:,5)) >= upperdepth); point_temp=nanmean(castprops(i_inrange,6)); point_sal=nanmean(castprops(i_inrange,7)); point_c=nanmean(castprops(i_inrange,8)); point_st=nanmean(castprops(i_inrange,9)); numdate=datenum(castprops(1,3),castprops(1,4),15); regcast=[castprops(1,1) castprops(1,2) castprops(1,3) castprops(1,4) mean([-abs(lowerdepth) -abs(upperdepth)]) point_temp point_sal point_c point_st]; allcasts=[allcasts; regcast]; end; placeholder=rownum; end; end; disp(['...',num2str(size(allcasts,1)),' found.']); if(qc) % Quality control phase f_qc=figure; tsdiagrm([min(allcasts(:,7))-1 max(allcasts(:,7))+1],[min(allcasts(:,6))-1 max(allcasts(:,6))+1],abs(allcasts(1,5))); hold on; plot(allcasts(:,7),allcasts(:,6),'b.'); disp(' '); disp('You have chosen to quality control the data. The T-S diagram shows temperature and salinity of all of the extracted data points--note that unlike a traditional T-S diagram, this plots shows only points from a single depth. You can eliminate outliers by choosing the upper left and lower right corners of a box enclosing the GOOD points. The mean T/S point is shown by the magenta *, and the 1x and 2x standard deviation are shown by the dashed magenta boxes. Watch the plot title for instructions.'); % Plot the mean temperature salinity point plot(nanmean(allcasts(:,7)), nanmean(allcasts(:,6)),'m*','markersize',10); % Plot temperature std lines plot([nanmean(allcasts(:,7))-nanstd(allcasts(:,7)) nanmean(allcasts(:,7))+nanstd(allcasts(:,7))],[nanmean(allcasts(:,6))+nanstd(allcasts(:,6)) nanmean(allcasts(:,6))+nanstd(allcasts(:,6))],'m--'); plot([nanmean(allcasts(:,7))-2*nanstd(allcasts(:,7)) nanmean(allcasts(:,7))+2*nanstd(allcasts(:,7))],[nanmean(allcasts(:,6))+2*nanstd(allcasts(:,6)) nanmean(allcasts(:,6))+2*nanstd(allcasts(:,6))],'m--'); plot([nanmean(allcasts(:,7))-nanstd(allcasts(:,7)) nanmean(allcasts(:,7))+nanstd(allcasts(:,7))],[nanmean(allcasts(:,6))-nanstd(allcasts(:,6)) nanmean(allcasts(:,6))-nanstd(allcasts(:,6))],'m--'); plot([nanmean(allcasts(:,7))-2*nanstd(allcasts(:,7)) nanmean(allcasts(:,7))+2*nanstd(allcasts(:,7))],[nanmean(allcasts(:,6))-2*nanstd(allcasts(:,6)) nanmean(allcasts(:,6))-2*nanstd(allcasts(:,6))],'m--'); % Plot salinity std lines plot([nanmean(allcasts(:,7))+nanstd(allcasts(:,7)) nanmean(allcasts(:,7))+nanstd(allcasts(:,7))],[nanmean(allcasts(:,6))-nanstd(allcasts(:,6)) nanmean(allcasts(:,6))+nanstd(allcasts(:,6))],'m--'); plot([nanmean(allcasts(:,7))+2*nanstd(allcasts(:,7)) nanmean(allcasts(:,7))+2*nanstd(allcasts(:,7))],[nanmean(allcasts(:,6))-2*nanstd(allcasts(:,6)) nanmean(allcasts(:,6))+2*nanstd(allcasts(:,6))],'m--'); plot([nanmean(allcasts(:,7))-nanstd(allcasts(:,7)) nanmean(allcasts(:,7))-nanstd(allcasts(:,7))],[nanmean(allcasts(:,6))-nanstd(allcasts(:,6)) nanmean(allcasts(:,6))+nanstd(allcasts(:,6))],'m--'); plot([nanmean(allcasts(:,7))-2*nanstd(allcasts(:,7)) nanmean(allcasts(:,7))-2*nanstd(allcasts(:,7))],[nanmean(allcasts(:,6))-2*nanstd(allcasts(:,6)) nanmean(allcasts(:,6))+2*nanstd(allcasts(:,6))],'m--'); title('Choose UPPER LEFT corner of box around points you wish to INCLUDE'); [s_ll,t_ul]=ginput(1); plot([s_ll s_ll],[min(allcasts(:,6))-1 t_ul],'g-'); plot([s_ll max(allcasts(:,7))+1],[t_ul t_ul],'g-'); title('Choose LOWER RIGHT corner of box around points you wish to INCLUDE'); [s_ul,t_ll]=ginput(1); % Find "good" and out-of-range points... i_good=find(allcasts(:,7) >= s_ll & allcasts(:,7) <= s_ul & allcasts(:,6) >= t_ll & allcasts(:,6) <= t_ul); i_bad=find(allcasts(:,7) < s_ll | allcasts(:,7) > s_ul | allcasts(:,6) < t_ll | allcasts(:,6) > t_ul); GOODPOINTS = allcasts(i_good,:); % Make new plot showing the eliminated points and new means/stds clf; tsdiagrm([min(allcasts(:,7))-1 max(allcasts(:,7))+1],[min(allcasts(:,6))-1 max(allcasts(:,6))+1],abs(allcasts(1,5))); hold on; plot(GOODPOINTS(:,7),GOODPOINTS(:,6),'b.'); % Plot the NEW mean temperature salinity point plot(nanmean(GOODPOINTS(:,7)), nanmean(GOODPOINTS(:,6)),'m*','markersize',10); % Plot NEW temperature std lines plot([nanmean(GOODPOINTS(:,7))-nanstd(GOODPOINTS(:,7)) nanmean(GOODPOINTS(:,7))+nanstd(GOODPOINTS(:,7))],[nanmean(GOODPOINTS(:,6))+nanstd(GOODPOINTS(:,6)) nanmean(GOODPOINTS(:,6))+nanstd(GOODPOINTS(:,6))],'m--'); plot([nanmean(GOODPOINTS(:,7))-2*nanstd(GOODPOINTS(:,7)) nanmean(GOODPOINTS(:,7))+2*nanstd(GOODPOINTS(:,7))],[nanmean(GOODPOINTS(:,6))+2*nanstd(GOODPOINTS(:,6)) nanmean(GOODPOINTS(:,6))+2*nanstd(GOODPOINTS(:,6))],'m--'); plot([nanmean(GOODPOINTS(:,7))-nanstd(GOODPOINTS(:,7)) nanmean(GOODPOINTS(:,7))+nanstd(GOODPOINTS(:,7))],[nanmean(GOODPOINTS(:,6))-nanstd(GOODPOINTS(:,6)) nanmean(GOODPOINTS(:,6))-nanstd(GOODPOINTS(:,6))],'m--'); plot([nanmean(GOODPOINTS(:,7))-2*nanstd(GOODPOINTS(:,7)) nanmean(GOODPOINTS(:,7))+2*nanstd(GOODPOINTS(:,7))],[nanmean(GOODPOINTS(:,6))-2*nanstd(GOODPOINTS(:,6)) nanmean(GOODPOINTS(:,6))-2*nanstd(GOODPOINTS(:,6))],'m--'); % Plot NEW salinity std lines plot([nanmean(GOODPOINTS(:,7))+nanstd(GOODPOINTS(:,7)) nanmean(GOODPOINTS(:,7))+nanstd(GOODPOINTS(:,7))],[nanmean(GOODPOINTS(:,6))-nanstd(GOODPOINTS(:,6)) nanmean(GOODPOINTS(:,6))+nanstd(GOODPOINTS(:,6))],'m--'); plot([nanmean(GOODPOINTS(:,7))+2*nanstd(GOODPOINTS(:,7)) nanmean(GOODPOINTS(:,7))+2*nanstd(GOODPOINTS(:,7))],[nanmean(GOODPOINTS(:,6))-2*nanstd(GOODPOINTS(:,6)) nanmean(GOODPOINTS(:,6))+2*nanstd(GOODPOINTS(:,6))],'m--'); plot([nanmean(GOODPOINTS(:,7))-nanstd(GOODPOINTS(:,7)) nanmean(GOODPOINTS(:,7))-nanstd(GOODPOINTS(:,7))],[nanmean(GOODPOINTS(:,6))-nanstd(GOODPOINTS(:,6)) nanmean(GOODPOINTS(:,6))+nanstd(GOODPOINTS(:,6))],'m--'); plot([nanmean(GOODPOINTS(:,7))-2*nanstd(GOODPOINTS(:,7)) nanmean(GOODPOINTS(:,7))-2*nanstd(GOODPOINTS(:,7))],[nanmean(GOODPOINTS(:,6))-2*nanstd(GOODPOINTS(:,6)) nanmean(GOODPOINTS(:,6))+2*nanstd(GOODPOINTS(:,6))],'m--'); plot([s_ll s_ll s_ul s_ul s_ll],[t_ll t_ul t_ul t_ll t_ll],'g-'); plot(allcasts(i_bad,7),allcasts(i_bad,6),'r.'); title('Blue points kept, red points eliminated'); disp(' '); disp(['Number of casts eliminated using T-S diagram: ',num2str(size(i_bad,1))]); disp(['Number eliminated because casts too shallow for chosen depth range: ',num2str(size(allcasts,1)-size(i_bad,1)-size(i_good,1))]); else % Eliminate Temperatures of NaN igoodtemps=find(isnan(allcasts(:,6))==0); GOOD3POINTS=allcasts(igoodtemps,:); % Eliminate Salinities of NaN igoodsals=find(isnan(GOOD3POINTS(:,7))==0); GOODPOINTS=GOOD3POINTS(igoodsals,:); disp(' '); disp(['Number eliminated because casts too shallow for chosen depth range: ',num2str(size(allcasts,1)-size(GOODPOINTS,1))]); end; % qc disp(' '); disp(['Total number of casts used: ',num2str(size(GOODPOINTS,1))]); if(plotyes>0) disp(' '); disp('You will now be shown a plot of your chosen averaging scheme.'); disp('The green points are the actual data points which will be averaged.'); disp('(i.e. stations that are too shallow/deep are NOT plotted.'); if radiusinc==0 tstr=['months: ',season,'; dep: ',num2str(upperdepth),'-',num2str(lowerdepth),'m; ','grid: ',num2str(gridspacing),'^{o}; radius: ',num2str(Ro),'km']; else tstr=['months: ',season,'; dep: ',num2str(upperdepth),'-',num2str(lowerdepth),'m; ','grid: ',num2str(gridspacing),'^{o}; radius: ',num2str(Ro),'km+']; end; f_grid=figure; % Grid setup m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); hold on; m_plot(GOODPOINTS(:,1),GOODPOINTS(:,2),'g.'); [cb,hb]=m_tbase('contour',[-cutofflowdepth -cutoffhighdepth],'k'); hold on; clabel(cb,hb,'color','black','fontsize',8); m_plot(X(:,5),Y(:,5),'k.'); hold on; m_plot(X(5,:),Y(5,:),'k.'); hold on; m_range_ring(X(:,5),Y(:,5),Ro,'color','red'); m_range_ring(X(5,:),Y(5,:),Ro,'color','red'); m_grid('box','fancy'); title(['Green=data; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'grid.png']); end; end; %%%%%%%%%%%%%%%%%%%%%%%%%%%% % Computing Planview means % %%%%%%%%%%%%%%%%%%%%%%%%%%%% disp(' '); disp('Computing planview means and standard deviations...'); NUMPOINTS=[]; RADIUS=[]; MEANTEMP=[]; MEANSAL=[]; MEANC=[]; MEANST=[]; STDTEMP=[]; STDSAL=[]; STDC=[]; STDST=[]; % Loop through every output grid node for loop1=1:size(newlat,2) for loop2=1:size(newlon,2) outputlat=Y(loop1,loop2); outputlon=X(loop1,loop2); % Compute the distance of every cast from this point RNG=deg2km(distance(ones(size(GOODPOINTS(:,2)))*outputlat,ones(size(GOODPOINTS(:,1)))*outputlon,GOODPOINTS(:,2),GOODPOINTS(:,1))); % Collect the ones within range (specified by Ro) Ronew=Ro; iin=find(RNG<=Ronew); if radiusinc > 0 NUMPOINTS(loop1,loop2)=size(iin,1); while length(iin) < minnumpts Ronew = Ronew + radiusinc; iin=find(RNG<=Ronew); end; RADIUS(loop1,loop2)=Ronew; f=0.54 + 0.46 * cos((RNG(iin)*pi)/Ronew); MEANTEMP(loop1,loop2)=nansum(f .* GOODPOINTS(iin,6))/nansum(f); junk=f .* (GOODPOINTS(iin,6)-MEANTEMP(loop1,loop2)); STDTEMP(loop1,loop2)=sqrt(nansum(junk.^2)/nansum(f)); MEANSAL(loop1,loop2)=nansum(f .* GOODPOINTS(iin,7))/nansum(f); junk=f .* (GOODPOINTS(iin,7)-MEANSAL(loop1,loop2)); STDSAL(loop1,loop2)=sqrt(nansum(junk.^2)/nansum(f)); MEANC(loop1,loop2)=nansum(f .* GOODPOINTS(iin,8))/nansum(f); junk=f .* (GOODPOINTS(iin,8)-MEANC(loop1,loop2)); STDC(loop1,loop2)=sqrt(nansum(junk.^2)/nansum(f)); MEANST(loop1,loop2)=nansum(f .* GOODPOINTS(iin,9))/nansum(f); junk=f .* (GOODPOINTS(iin,9)-MEANST(loop1,loop2)); STDST(loop1,loop2)=sqrt(nansum(junk.^2)/nansum(f)); elseif length(iin) < minnumpts NUMPOINTS(loop1,loop2)=NaN; RADIUS(loop1,loop2)=Ronew; MEANTEMP(loop1,loop2)=NaN; STDTEMP(loop1,loop2)=NaN; MEANSAL(loop1,loop2)=NaN; STDSAL(loop1,loop2)=NaN; MEANC(loop1,loop2)=NaN; STDC(loop1,loop2)=NaN; MEANST(loop1,loop2)=NaN; STDST(loop1,loop2)=NaN; else NUMPOINTS(loop1,loop2)=size(iin,1); f=0.54 + 0.46 * cos((RNG(iin)*pi)/Ronew); MEANTEMP(loop1,loop2)=nansum(f .* GOODPOINTS(iin,6))/nansum(f); junk=f .* (GOODPOINTS(iin,6)-MEANTEMP(loop1,loop2)); STDTEMP(loop1,loop2)=sqrt(nansum(junk.^2)/nansum(f)); MEANSAL(loop1,loop2)=nansum(f .* GOODPOINTS(iin,7))/nansum(f); junk=f .* (GOODPOINTS(iin,7)-MEANSAL(loop1,loop2)); STDSAL(loop1,loop2)=sqrt(nansum(junk.^2)/nansum(f)); MEANC(loop1,loop2)=nansum(f .* GOODPOINTS(iin,8))/nansum(f); junk=f .* (GOODPOINTS(iin,8)-MEANC(loop1,loop2)); STDC(loop1,loop2)=sqrt(nansum(junk.^2)/nansum(f)); MEANST(loop1,loop2)=nansum(f .* GOODPOINTS(iin,9))/nansum(f); junk=f .* (GOODPOINTS(iin,9)-MEANST(loop1,loop2)); STDST(loop1,loop2)=sqrt(nansum(junk.^2)/nansum(f)); end; end; end; % Need to NaN out land pixels outputnodeelevs=interp2(LON,LAT,data,X,Y); ibad=find(outputnodeelevs>-1*abs(upperdepth)); NUMPOINTS(ibad)=NaN; MEANTEMP(ibad)=NaN; MEANSAL(ibad)=NaN; MEANC(ibad)=NaN; MEANST(ibad)=NaN; STDTEMP(ibad)=NaN; STDSAL(ibad)=NaN; STDC(ibad)=NaN; STDST(ibad)=NaN; if(plotyes>0) disp(' '); disp('Plotting...'); if radiusinc > 0 RADIUS(ibad)=NaN; f0=figure; % Final search radius m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,RADIUS,16); caxis([Ro max(max(RADIUS))]); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); title(['RADIUS; ',tstr]); set(findobj('tag','m_grid_color'),'facecolor','none'); if plotyes==2 print('-dpng',[saveprefix 'radius.png']); end; end; f1=figure; % Number of points per output grid node m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,NUMPOINTS,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['NUMPOINTS; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'numpoints.png']); end; f2=figure; % Mean temperature m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,MEANTEMP,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['MEAN TEMP; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'meantemp.png']); end; f3=figure; % Standard deviation temperature m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,STDTEMP,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['STD TEMP; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'stdtemp.png']); end; f4=figure; % Mean salinity m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,MEANSAL,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['MEAN SAL; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'meansal.png']); end; f5=figure; % Standard deviation salinity m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,STDSAL,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['STD SAL; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'stdsal.png']); end; f6=figure; % Mean Sigma-t m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,MEANST-1000,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['MEAN SIGMAT; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'meanst.png']); end; f7=figure; % Standard deviation sigma-t m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,STDST,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['STD SIGMAT; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'stdst.png']); end; f8=figure; % Mean sound speed m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,MEANC,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['MEAN C; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'meanc.png']); end; f9=figure; % Standard deviation sound speed m_proj('mercator','lon',[westbound eastbound],'lat',[southbound northbound]); hold on; [c,h]=m_contourf(X,Y,STDC,16); colorbar; for fixloop=1:length(h) set(h(fixloop),'linestyle','none'); end; m_gshhs_i('patch',[.7 .7 .7]); m_grid('box','fancy'); set(findobj('tag','m_grid_color'),'facecolor','none'); title(['STD C; ',tstr]); if plotyes==2 print('-dpng',[saveprefix 'stdc.png']); end; end; % plotyes disp(' '); disp('Program finished.'); disp(' ');