/[MITgcm]/MITgcm_contrib/eh3/regrid/scrip_examples/remap_grid_CS32.m
ViewVC logotype

Annotation of /MITgcm_contrib/eh3/regrid/scrip_examples/remap_grid_CS32.m

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


Revision 1.1 - (hide annotations) (download)
Thu Aug 10 05:00:16 2006 UTC (19 years ago) by edhill
Branch: MAIN
CVS Tags: HEAD
initial ci

1 edhill 1.1 %
2     % %Id: $
3     %
4     % Ed Hill
5     %
6     % Create a netCDF file defining a CS32 grid for SCRIP.
7    
8     clear all
9     close all
10    
11     debug_lev = 0;
12    
13     %==================================================================
14     % Read the tile00?.mitgrid files
15     gvars = { 'XC','YC','DXF','DYF','RA','XG','YG','DXV', ...
16     'DYU','RAZ','DXC','DYC','RAW','RAS','DXG','DYG' };
17     ne = 32;
18     nep1 = ne + 1;
19     iface = 1;
20     for iface = 1:6
21     fname = sprintf('CS32/tile%03d.mitgrid', iface);
22     gid = fopen(fname, 'r', 'ieee-be');
23     tmp = reshape(fread(gid,inf,'real*8',0,'ieee-be'),[nep1,nep1,16]);
24     fclose(gid);
25     % surf(tmp(:,:,1)), view(2), shading interp
26     % for jj = 1:length(gvars)
27     for jj = 1:7
28     comm = sprintf('%s(:,:,%d) = tmp(:,:,%d);', ...
29     [gvars{jj}], iface, jj);
30     eval(comm);
31     end
32     end
33     % surf(XC(:,:,1)), view(2), shading interp
34     % subplot(2,1,1), a = [1:10]; surf(XC(a,a,1)), view(2)
35     % subplot(2,1,2), a = [(nep1-10):nep1]; surf(XC(a,a,1)), view(2)
36     % surf(YC(:,:,1)), view(2), shading interp
37     % surf(XG(:,:,1)), view(2), shading interp
38     % surf(YG(:,:,1)), view(2), shading interp
39     %==================================================================
40    
41    
42     grid_size = ne * ne * 6;
43     grid_rank = 2;
44     grid_corners = 4;
45    
46     ! rm -f remap_grid_CS32.cdf
47    
48     nc = netcdf(['remap_grid_CS32.cdf'], 'clobber');
49     nc.title = 'Global CS32 grid';
50     nc.author = 'Ed Hill';
51     nc('grid_size') = grid_size;
52     nc('grid_rank') = grid_rank;
53     nc('grid_corners') = grid_corners;
54    
55     nc{'grid_dims'} = ncint('grid_rank');
56     nc{'grid_dims'}(:) = [ ne 6*ne ];
57    
58     nc{'grid_imask'} = ncint('grid_size');
59     nc{'grid_imask'}(:) = ones(length(grid_size),1);
60    
61     nc{'grid_center_lat'} = ncdouble('grid_size');
62     nc{'grid_center_lon'} = ncdouble('grid_size');
63     nc{'grid_center_lat'}.units = 'degrees';
64     nc{'grid_center_lon'}.units = 'degrees';
65     nc{'grid_corner_lat'} = ncdouble('grid_size', 'grid_corners');
66     nc{'grid_corner_lon'} = ncdouble('grid_size', 'grid_corners');
67     nc{'grid_corner_lat'}.units = 'degrees';
68     nc{'grid_corner_lon'}.units = 'degrees';
69    
70     cen_lat = zeros(grid_size,1);
71     cen_lon = zeros(grid_size,1);
72     cor_lat = zeros(grid_size,4);
73     cor_lon = zeros(grid_size,4);
74     j = 0;
75     for k = 1:6
76     for iy = 1:ne
77     for ix = 1:ne
78     j = j + 1;
79     cen_lat(j) = YC(ix,iy,k);
80     cen_lon(j) = XC(ix,iy,k);
81     % Note: grid points *MUST* be ordered CCW and here
82     % we start from the lower-left corner
83     a = [ YG(ix,iy,k) YG(ix+1,iy,k) YG(ix+1,iy+1,k) YG(ix,iy+1,k) ];
84     b = [ XG(ix,iy,k) XG(ix+1,iy,k) XG(ix+1,iy+1,k) XG(ix,iy+1,k) ];
85     % fixup the poles
86     for ii = 1:4
87     if (abs(b(ii)-XC(ix,iy,k)) > 150)
88     b(ii) = b(ii) + sign(XC(ix,iy,k))*360;
89     end
90     if (90 - abs(a(ii))) < 0.01
91     b(ii) = XC(ix,iy,k);
92     end
93     end
94     cor_lat(j,:) = a;
95     cor_lon(j,:) = b;
96     % if debug_lev > 1
97     % pause(1)
98     % hold on
99     % plot([cor_lon(j,:) cor_lon(j,1)],...
100     % [cor_lat(j,:) cor_lat(j,1)],'o-b')
101     % plot(cen_lon(j),cen_lat(j),'xr')
102     % hold off
103     % end
104     end
105     end
106     end
107     nc{'grid_center_lat'}(:) = cen_lat;
108     nc{'grid_center_lon'}(:) = cen_lon;
109     nc{'grid_corner_lat'}(:) = cor_lat;
110     nc{'grid_corner_lon'}(:) = cor_lon;
111    
112     close(nc)
113    

  ViewVC Help
Powered by ViewVC 1.1.22