Home > atmlab > arts > various > arts_plot_atmgrids.m

arts_plot_atmgrids

PURPOSE ^

ARTS_PLOT_ATMGRIDS Makes figure for 1D-3D with atmospheric grids.

SYNOPSIS ^

function [h,hs] = arts_plot_atmgrids(dim,lat_grid,lon_grid,z_field,r_geoid,z_ground,cb_lims,mark_fields)

DESCRIPTION ^

 ARTS_PLOT_ATMGRIDS   Makes figure for 1D-3D with atmospheric grids.

     The function produces figures to visualize the atmsopheric grids,
     the geoid, the ground, the cloud box and the atmospheric fields.

     The variables work as their ARTS counterparts. The exception is for
     1D where *lat_grid* shall contain two values, which determines how
     long (in latitude) the plotted atmosphere shall be extended.

     Examples on how to use this function is found in ARTS, at
     doc/uguide/Figs/fm_definitions/mkfigs_atm_dims.m

 FORMAT   [h,hs] = arts_plot_atmgrids(dim,lat_grid,lon_grid,z_field,
                                       r_geoid,[z_ground,cb_lims,mark_fields])
        
 OUT   h         Handles to plot symbols.
       hs        Text suitable to put in legend.
 IN    dim       Atmospheric dimension (1-3).
       lat_grid  As the WSV with the same name, except for 1D. See above.
       lon_grid  As the WSV with the same name.
       z_field   As the WSV with the same name.
       r_geoid   As the WSV with the same name.
 OPT   z_ground  As the WSV with the same name. If the variable is empty,
                 which is default, the geoid and the ground are not plotted.
       cb_lims   As the WSV cloudbox_limits. The cloud box is turned off
                 with cb_lims=[], which is also default.
       mark_fields  Flag to indicate atmospheric fields with black dot at
                 crossings of grids. Default is 0.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

DOWNLOAD ^

arts_plot_atmgrids.m

SOURCE CODE ^

0001 % ARTS_PLOT_ATMGRIDS   Makes figure for 1D-3D with atmospheric grids.
0002 %
0003 %     The function produces figures to visualize the atmsopheric grids,
0004 %     the geoid, the ground, the cloud box and the atmospheric fields.
0005 %
0006 %     The variables work as their ARTS counterparts. The exception is for
0007 %     1D where *lat_grid* shall contain two values, which determines how
0008 %     long (in latitude) the plotted atmosphere shall be extended.
0009 %
0010 %     Examples on how to use this function is found in ARTS, at
0011 %     doc/uguide/Figs/fm_definitions/mkfigs_atm_dims.m
0012 %
0013 % FORMAT   [h,hs] = arts_plot_atmgrids(dim,lat_grid,lon_grid,z_field,
0014 %                                       r_geoid,[z_ground,cb_lims,mark_fields])
0015 %
0016 % OUT   h         Handles to plot symbols.
0017 %       hs        Text suitable to put in legend.
0018 % IN    dim       Atmospheric dimension (1-3).
0019 %       lat_grid  As the WSV with the same name, except for 1D. See above.
0020 %       lon_grid  As the WSV with the same name.
0021 %       z_field   As the WSV with the same name.
0022 %       r_geoid   As the WSV with the same name.
0023 % OPT   z_ground  As the WSV with the same name. If the variable is empty,
0024 %                 which is default, the geoid and the ground are not plotted.
0025 %       cb_lims   As the WSV cloudbox_limits. The cloud box is turned off
0026 %                 with cb_lims=[], which is also default.
0027 %       mark_fields  Flag to indicate atmospheric fields with black dot at
0028 %                 crossings of grids. Default is 0.
0029 
0030 % 2003-03-04   Created by Patrick Eriksson.
0031 
0032 function [h,hs] = arts_plot_atmgrids(dim,lat_grid,lon_grid,z_field,...
0033                                           r_geoid,z_ground,cb_lims,mark_fields)
0034 
0035 %= Check input                                                   %&%
0036 %                                                                %&%
0037 rqre_nargin(5,nargin);                                           %&%
0038 %                                                                %&%
0039 if ~isscalar( dim )  |  dim < 1  |  dim > 3                      %&%
0040   error( ['The atmospheric dimensionality must be a ',...        %&%
0041                                   'scalar between 1 and 3.' ] ); %&%
0042 end                                                              %&%
0043 %                                                                %&%
0044 if dim == 1  & ( ~isvector( lat_grid )   |  ...                  %&%
0045                                        length( lat_grid ) ~= 2 ) %&%
0046   error('For 1D, *lat_gid* must be a vector of length 2.');      %&%
0047 end                                                              %&%
0048                                                                  %&%
0049 
0050 
0051 if nargin < 6
0052   z_ground = [];
0053 end
0054 %
0055 if nargin < 7
0056   cb_lims = [];
0057 end
0058 %
0059 if nargin < 8
0060   mark_fields = 0;
0061 end
0062 
0063 
0064 
0065 %= Check if plot is hold
0066 %
0067 figure(gcf);
0068 figstat = ishold;
0069 %
0070 if ~ishold
0071   clf
0072 end
0073 
0074 
0075 %= Plotting symbols
0076 %
0077 ps_grids     = 'k:';
0078 ps_ground{1} = 'r-';
0079 ps_ground{2} = 'LineWidth';
0080 ps_ground{3} = 2;
0081 ps_geoid     = 'b-.';
0082 ps_cloudb    = 0.7*ones(1,3);
0083 ps_field     = 'k.';
0084 %
0085 hs{1} = 'Atmospheric grids';
0086 hs{2} = 'Ground';
0087 hs{3} = 'Geoid';
0088 if ~isempty( cb_lims )
0089   hs{4} = 'Cloud box';
0090   do_cb = 1;
0091 else
0092   do_cb = 0;
0093 end
0094 if mark_fields
0095   hs{4+do_cb} = 'Atmospheric field';
0096 end
0097 
0098 
0099 nz = size(z_field,1);
0100 
0101 
0102 switch dim
0103 
0104   %- 1D
0105   case 1
0106     %
0107     if ~isempty( cb_lims )
0108       [x1,y1] = surf_segment2D( lat_grid, ... 
0109                                         [r_geoid(1)+z_field(cb_lims{1}+1), ...
0110                                          r_geoid(1)+z_field(cb_lims{1}+1) ] );
0111       [x2,y2] = surf_segment2D( lat_grid, ...
0112                                         [r_geoid(1)+z_field(cb_lims{2}+1), ...
0113                                          r_geoid(1)+z_field(cb_lims{2}+1) ] );
0114       h(4) = patch( [x1 fliplr(x2)], [y1 fliplr(y2)], ps_cloudb);
0115       set( h(4), 'EdgeColor', ps_cloudb );
0116       hold on
0117     end
0118     %
0119     for j = 1 : nz
0120       [x,y] = surf_segment2D( lat_grid, ...
0121                             [ r_geoid(1)+z_field(j), r_geoid(1)+z_field(j) ] );
0122       h(1) = plot( x, y, ps_grids );
0123       hold on
0124       %
0125       if mark_fields
0126         i0 = round( length(x)/2 );
0127         h(4+do_cb) = plot( x(i0), y(i0), ps_field );
0128       end
0129 
0130     end
0131     %
0132     if ~isempty( z_ground )
0133       [x,y] = surf_segment2D( lat_grid, ...
0134                           [ r_geoid(1)+z_ground(1), r_geoid(1)+z_ground(1) ] );
0135       h(2) = plot( x, y, ps_ground{:} );
0136       %
0137       [x,y] = surf_segment2D( lat_grid, [ r_geoid(1), r_geoid(1) ] );
0138       h(3) = plot( x, y, ps_geoid );
0139     end
0140 
0141   %- 2D
0142   case 2
0143     %
0144     if ~isempty( cb_lims )
0145       ind     = (cb_lims{3}+1):(cb_lims{4}+1);
0146       [x1,y1] = surf_segment2D( lat_grid(ind), ...
0147                                      r_geoid(ind)+z_field(cb_lims{1}+1,ind)' );
0148       [x2,y2] = surf_segment2D( lat_grid(ind), ...
0149                                      r_geoid(ind)+z_field(cb_lims{2}+1,ind)' );
0150       h(4) = patch( [x1 fliplr(x2)], [y1 fliplr(y2)], ps_cloudb );
0151       set( h(4), 'EdgeColor', ps_cloudb );
0152       hold on
0153     end
0154     %
0155     for j = 1 : nz
0156       [x,y] = surf_segment2D( lat_grid, r_geoid+z_field(j,:)' );
0157       h(1) = plot( x, y, ps_grids );
0158       hold on
0159     end
0160     %
0161     for j = 1 : length( lat_grid )
0162       [x,y] = atmplot_pol2cart( [ r_geoid(j)+z_field(1,j), ...
0163                            r_geoid(j)+z_field(nz,j) ], lat_grid(j)*ones(1,2) );
0164       plot( x, y, ps_grids );
0165       %
0166       if mark_fields
0167         for k = 1 : nz
0168           [x,y] = atmplot_pol2cart( r_geoid(j)+z_field(k,j), lat_grid(j) );
0169           h(4+do_cb) = plot( x, y, ps_field );
0170         end
0171       end
0172     end
0173     %
0174     if ~isempty( z_ground )
0175       [x,y] = surf_segment2D( lat_grid, r_geoid+z_ground );
0176       h(2) = plot( x, y, ps_ground{:} );
0177       %
0178       [x,y] = surf_segment2D( lat_grid, r_geoid );
0179       h(3) = plot( x, y, ps_geoid );
0180     end
0181 
0182   %- 3D
0183   case 3
0184     %
0185     for j = 1 : length( lat_grid )
0186       for k = 1 : length( lon_grid )
0187         [x,y,z] = atmplot_sph2cart( [ r_geoid(j,k)+z_field(1,j,k), ...
0188                                 r_geoid(j,k)+z_field(nz,j,k) ], ...
0189                     lat_grid(j)*ones(1,2), lon_grid(k)*ones(1,2) );
0190           h(1) = plot3( z, x, y, ps_grids );
0191         hold on
0192           %
0193           if mark_fields
0194             for l = 1 : nz
0195               [x,y,z] = atmplot_sph2cart( r_geoid(j,k)+z_field(l,j,k), ...
0196                                                     lat_grid(j), lon_grid(k) );
0197               h(4+do_cb) = plot3( z, x, y, ps_field );
0198             end
0199           end
0200       end
0201     end
0202     %
0203     for j = 1 : length( lat_grid )
0204       if ~isempty( z_ground )
0205         [x,y,z] = surf_segment3D( lat_grid(j), lon_grid, ...
0206                                                r_geoid(j,:)+z_ground(j,:,:) );
0207         h(2) = plot3( z, x, y, ps_ground{:} );
0208         %
0209         [x,y,z] = surf_segment3D( lat_grid(j), lon_grid, r_geoid(j,:) );
0210         h(3) = plot3( z, x, y, ps_geoid );
0211       end
0212       %
0213       for k = 1 : nz
0214         [x,y,z] = surf_segment3D( lat_grid(j), lon_grid, ...
0215                                       r_geoid(j,:)+squeeze(z_field(k,j,:))' );
0216         plot3( z, x, y, ps_grids );
0217       end
0218     end
0219     %
0220     for j = 1 : length( lon_grid )
0221       if ~isempty( z_ground )
0222 
0223         [x,y,z] = surf_segment3D( lat_grid, lon_grid(j), ...
0224                                                   r_geoid(:,j)+z_ground(:,j) );
0225         plot3( z, x, y, ps_ground{:} );
0226         %
0227         [x,y,z] = surf_segment3D( lat_grid, lon_grid(j), r_geoid(:,j) );
0228         plot3( z, x, y, ps_geoid );
0229       end
0230       %
0231       for k = 1 : nz
0232         [x,y,z] = surf_segment3D( lat_grid, lon_grid(j), ...
0233                                        r_geoid(:,j)+squeeze(z_field(k,:,j))' );
0234         plot3( z, x, y, ps_grids );
0235       end
0236     end
0237     %
0238     set( gca, 'Ydir', 'rev' );
0239     view([37 30]);
0240 
0241 end
0242 
0243 
0244 %= Set equal axis
0245 %
0246 axis equal
0247 
0248 
0249 %= Set hold off if this was status at the start
0250 %
0251 if ~figstat
0252   hold off
0253 end
0254 
0255 
0256 
0257 function [x,y] = surf_segment2D(lat,r)
0258 
0259   % Nominal latitude step length
0260   nom_dlat = 0.2;
0261 
0262   x = [];
0263   y = [];
0264 
0265   for j = 2 : length( lat )
0266 
0267     lat1 = lat(j-1);
0268     lat3 = lat(j);
0269 
0270     lats = linspace( lat1, lat3, max([ 2, ceil( (lat3-lat1) / nom_dlat ) ] ) );
0271 
0272     slope = ( r(j) - r(j-1) ) / ( lat3 - lat1 );
0273 
0274     rs    = r(j-1) + slope * ( lats - lat1 );
0275       
0276     [xv,yv] = atmplot_pol2cart( rs, lats );
0277 
0278     x = [ x xv ];
0279     y = [ y yv ];
0280 
0281   end
0282 return
0283 
0284 
0285 
0286 function [x,y,z] = surf_segment3D(lat,lon,r)
0287 
0288   % Nominal latitude step length
0289   nom_dang = 0.5;
0290 
0291   x = [];
0292   y = [];
0293   z = [];
0294 
0295   if length(lon) == 1
0296     for j = 2 : length( lat )
0297       lat1 = lat(j-1);
0298       lat3 = lat(j);
0299       lats = linspace( lat1, lat3, max([ 2, ceil( (lat3-lat1) / nom_dang ) ]));
0300       slope = ( r(j) - r(j-1) ) / ( lat3 - lat1 );
0301       rs    = r(j-1) + slope * ( lats - lat1 );
0302       [xv,yv,zv] = atmplot_sph2cart( rs, lats, lon );
0303       x = [ x xv ];
0304       y = [ y yv ];
0305       z = [ z zv ];
0306     end
0307   else
0308     for j = 2 : length( lon )
0309       lon1 = lon(j-1);
0310       lon3 = lon(j);
0311       lons = linspace( lon1, lon3, max([ 2, ceil( (lon3-lon1) / nom_dang ) ]));
0312       slope = ( r(j) - r(j-1) ) / ( lon3 - lon1 );
0313       rs    = r(j-1) + slope * ( lons - lon1 );
0314       [xv,yv,zv] = atmplot_sph2cart( rs, lat, lons );
0315       x = [ x xv ];
0316       y = [ y yv ];
0317       z = [ z zv ];
0318     end
0319   end
0320 return

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005