9 integer*4,
parameter :: idim=43200
10 integer*4,
parameter :: jdim=21600
11 integer*4,
parameter :: idim_p1=43201
12 integer*4,
parameter :: jdim_p1=21601
14 character(len=150) :: filenetcdf, fileraw
16 integer :: i, istat, ncid, status, dim_i, dim_j
17 integer :: dim_ip1, dim_jp1
18 integer :: id_lon, id_lat, id_data
19 integer :: id_lat_corner, id_lon_corner
21 integer(kind=2),
allocatable :: topo(:,:)
23 real(kind=8),
allocatable :: lats(:), lons(:)
24 real(kind=8),
allocatable :: lats_corner(:), lons_corner(:)
25 real(kind=8) :: lat11, lon11, dx, dy
30 lat11 = 90.0_8 + dy*0.5_8
31 lon11 = -180.0_8 + dx*0.5_8
33 allocate(lons(idim),lats(jdim),topo(idim,jdim))
34 allocate(lons_corner(idim_p1),lats_corner(jdim_p1))
37 lons(i) =
real((i-1),8) * dx + lon11
38 print*,
'lon ',i,lons(i)
42 lats(i) =
real((i-1),8) * dy + lat11
43 print*,
'lat ',i,lats(i)
50 lons_corner(i) =
real((i-1),8) * dx + lon11
51 print*,
'lon_corner ',i,lons_corner(i)
55 lats_corner(i) =
real((i-1),8) * dy + lat11
56 print*,
'lat_corner ',i,lats_corner(i)
59 fileraw=
"/scratch1/NCEPDEV/global/glopara/fix/raw/orog/gmted2010.30sec.int"
61 open(11, file=trim(fileraw), access=
'direct', recl=idim*jdim*2)
62 read(11, rec=1, iostat=istat) topo
63 if (istat /= 0) stop 99
66 print*,
'topo ', maxval(topo),minval(topo)
68 filenetcdf=
"./topography.gmted2010.30s.nc"
70 print*,
"- CREATE FILE: ", trim(filenetcdf)
71 status=nf90_create(filenetcdf, ior(nf90_netcdf4,nf90_classic_model), ncid)
72 if (status /= nf90_noerr) stop 1
74 status=nf90_def_dim(ncid,
'idim', idim, dim_i)
75 if (status /= nf90_noerr) stop 3
77 status=nf90_def_dim(ncid,
'jdim', jdim, dim_j)
78 if (status /= nf90_noerr) stop 2
80 status=nf90_def_dim(ncid,
'idim_p1', (idim+1), dim_ip1)
81 if (status /= nf90_noerr) stop 4
83 status=nf90_def_dim(ncid,
'jdim_p1', (jdim+1), dim_jp1)
84 if (status /= nf90_noerr) stop 5
86 status=nf90_put_att(ncid, nf90_global,
'source',
'USGS GMTED2010 TOPOGRAPHY DATA')
87 if (status /= nf90_noerr) stop 6
89 status=nf90_put_att(ncid, nf90_global,
'projection',
'regular lat/lon')
90 if (status /= nf90_noerr) stop 67
92 status=nf90_def_var(ncid,
'lat', nf90_double, dim_j, id_lat)
93 if (status /= nf90_noerr) stop 17
95 status=nf90_put_att(ncid, id_lat,
'long_name',
'grid cell center latitude')
96 if (status /= nf90_noerr) stop 10
98 status=nf90_put_att(ncid, id_lat,
'units',
'degrees')
99 if (status /= nf90_noerr) stop 85
101 status=nf90_def_var(ncid,
'lat_corner', nf90_double, dim_jp1, id_lat_corner)
102 if (status /= nf90_noerr) stop 37
104 status=nf90_put_att(ncid, id_lat_corner,
'long_name',
'grid cell corner latitude')
105 if (status /= nf90_noerr) stop 38
107 status=nf90_put_att(ncid, id_lat_corner,
'units',
'degrees')
108 if (status /= nf90_noerr) stop 86
110 status=nf90_def_var(ncid,
'lon', nf90_double, dim_i, id_lon)
111 if (status /= nf90_noerr) stop 16
113 status=nf90_put_att(ncid, id_lon,
'long_name',
'grid cell center longitude')
114 if (status /= nf90_noerr) stop 10
116 status=nf90_put_att(ncid, id_lon,
'units',
'degrees')
117 if (status /= nf90_noerr) stop 87
119 status=nf90_def_var(ncid,
'lon_corner', nf90_double, dim_ip1, id_lon_corner)
120 if (status /= nf90_noerr) stop 16
122 status=nf90_put_att(ncid, id_lon_corner,
'long_name',
'grid cell corner longitude')
123 if (status /= nf90_noerr) stop 40
125 status=nf90_put_att(ncid, id_lon_corner,
'units',
'degrees')
126 if (status /= nf90_noerr) stop 88
128 status=nf90_def_var(ncid,
'topo', nf90_short, (/dim_i,dim_j/), id_data)
129 if (status /= nf90_noerr) stop 20
131 status=nf90_put_att(ncid, id_data,
'long_name',
'topography')
132 if (status /= nf90_noerr) stop 65
134 status=nf90_put_att(ncid, id_data,
'units',
'meters')
135 if (status /= nf90_noerr) stop 55
137 status=nf90_enddef(ncid)
138 if (status /= nf90_noerr) stop 22
140 status=nf90_put_var(ncid, id_lon, lons)
141 if (status /= nf90_noerr) stop 19
143 status=nf90_put_var(ncid, id_lon_corner, lons_corner)
144 if (status /= nf90_noerr) stop 59
146 status=nf90_put_var(ncid, id_lat, lats)
147 if (status /= nf90_noerr) stop 20
149 status=nf90_put_var(ncid, id_lat_corner, lats_corner)
150 if (status /= nf90_noerr) stop 57
152 status=nf90_put_var(ncid, id_data, topo)
153 if (status /= nf90_noerr) stop 24
155 status=nf90_close(ncid)
157 end program topo_netcdf