IODA Bundle
gnssro_gsidiag2ioda.f90
Go to the documentation of this file.
1 !--- Hailing Zhang 2021-05-20
2 !--- This is to further process the GSI Diag netcdf file
3 ! from the JCSDA GSI branch feature/jedi_gdas
4 ! it 1) reorders the randomly distributed observations
5 ! and output observations in profile record numbers,
6 ! 2) adds variable descriptions, and
7 ! 3) output observation and geoval(optional) files seperately.
8 !-----------------------------------------------------------------
9 !-- TO COMPILE standalone executable
10 !ifort -o gnssro_gsidiag2ioda.exe
11 !gnssro_gsidiag2ioda.f90 -L${NETCDF}/lib -lnetcdf -lnetcdff
12 !-lm -I${NETCDF}/include -I${NETCDF}/include
13 !-- TO USE
14 ! gnssro_gsidiag2ioda.exe yyyymmddhh $gsidiag_input 1 !!!! 1
15 ! to write geoval/0 not
16 
17 program gnssro_gsidiag_full2small
18 use netcdf
19 implicit none
20 integer, parameter :: i_kind = selected_int_kind(8)
21 integer, parameter :: r_single = selected_real_kind(4)
22 
23 character(len=256) :: filename_in, obsname_out, geoname_out
24 integer :: ncid_in, ncid_out, ncid_out2
25 integer :: nobs_dimid, nlevs_dimid, nlevs1_dimid
26 integer :: nlevs, nlevs1, nobs, nlocs, nrecs
27 integer :: nlocs_dimid,nvars_dimid,nrecs_dimid, nlocs_dimid2
28 integer :: varid_meta,varid_geoval_nlevs,varid_geoval_nlevs1
29 integer :: varid_lat, varid_lon, varid_time
30 integer :: varid_lat2,varid_lon2,varid_time2
31 integer :: varid_geoid, varid_rfict
32 integer :: varid_rec, varid_said,varid_siid,varid_ogce,varid_ptid,varid_sclf,varid_asce
33 integer :: varid_imph,varid_impp, varid_azim
34 integer :: varid_bnd,varid_preqc, varid_gsihofx,varid_obserr,varid_obserr_adjust,varid_obserr_final
35 integer :: varid_geo_virtemp,varid_geo_airtemp,varid_geo_pres,varid_geo_shum,varid_geo_geop
36 integer :: varid_geo_pres1, varid_geo_geop1
37 integer :: varid_geo_geop_sfc
38 integer :: istart(2), icount(2)
39 real(r_single) :: zero_single
40 real(r_single), parameter :: huge_single = huge(zero_single) -1
41 character(len=10) :: anatime
42 character(len=1) :: geovalwrite
43 integer(i_kind) :: anatime_i
44 
45 type gpsro_type
46  integer, allocatable :: said(:)
47  integer, allocatable :: siid(:)
48  integer, allocatable :: ogce(:)
49  integer, allocatable :: sclf(:)
50  integer, allocatable :: ptid(:)
51  integer, allocatable :: asce(:)
52  integer, allocatable :: record_number(:)
53 
54  real(r_single), allocatable :: lat(:)
55  real(r_single), allocatable :: lon(:)
56  real(r_single), allocatable :: time(:)
57  real(r_single), allocatable :: geoid(:)
58  real(r_single), allocatable :: rfict(:)
59  real(r_single), allocatable :: azim(:)
60  real(r_single), allocatable :: impact_height(:)
61  real(r_single), allocatable :: impact_parameter(:)
62  real(r_single), allocatable :: bnd_obs(:)
63  integer, allocatable :: bnd_preqc(:)
64  real(r_single), allocatable :: bnd_obserr(:)
65  real(r_single), allocatable :: bnd_gsihofx(:)
66  real(r_single), allocatable :: bnd_obserr_final(:)
67  real(r_single), allocatable :: bnd_obserr_adjust(:)
68  real(r_single), allocatable :: surface_geopotential_height(:)
69  real(r_single), allocatable :: air_temperature(:,:)
70  real(r_single), allocatable :: virtual_temperature(:,:)
71  real(r_single), allocatable :: specific_humidity(:,:)
72  real(r_single), allocatable :: geopotential_height(:,:)
73  real(r_single), allocatable :: air_pressure(:,:)
74  real(r_single), allocatable :: geopotential_height_levels(:,:)
75  real(r_single), allocatable :: air_pressure_levels(:,:)
76 end type gpsro_type
77 
78 type(gpsro_type) :: gpsro_data
79 type(gpsro_type) :: gpsro_subset
80 type(gpsro_type) :: gpsro_reorder
81 integer :: nsub, i, j, k,substart, subend,narg
82 
83 narg=command_argument_count()
84 call getarg(1,anatime)
85 call getarg(2,filename_in)
86 ! optional 3rd argument: "1" means geoval will be written
87 if (narg == 3) then
88  call getarg(3,geovalwrite)
89 else
90  geovalwrite="1"
91 end if
92 read(anatime,'(i10)') anatime_i
93 
94 obsname_out = './gnssro_obs_'//trim(anatime)//'.nc4'
95 
96 call check(nf90_open(trim(filename_in), nf90_nowrite, ncid_in))
97 call check(nf90_inq_dimid(ncid_in, "nobs", nobs_dimid))
98 call check(nf90_inq_dimid(ncid_in, "air_temperature_arr_dim", nlevs_dimid))
99 call check(nf90_inq_dimid(ncid_in, "air_pressure_levels_arr_dim", nlevs1_dimid))
100 call check(nf90_inquire_dimension(ncid_in, nobs_dimid, len=nobs))
101 call check(nf90_inquire_dimension(ncid_in, nlevs_dimid, len=nlevs))
102 call check(nf90_inquire_dimension(ncid_in, nlevs1_dimid, len=nlevs1))
103 
104 allocate(gpsro_data%lat(nobs))
105 allocate(gpsro_data%lon(nobs))
106 allocate(gpsro_data%time(nobs))
107 allocate(gpsro_data%asce(nobs))
108 allocate(gpsro_data%azim(nobs))
109 allocate(gpsro_data%geoid(nobs))
110 allocate(gpsro_data%rfict(nobs))
111 allocate(gpsro_data%impact_height(nobs))
112 allocate(gpsro_data%impact_parameter(nobs))
113 allocate(gpsro_data%said(nobs))
114 allocate(gpsro_data%siid(nobs))
115 allocate(gpsro_data%ogce(nobs))
116 allocate(gpsro_data%ptid(nobs))
117 allocate(gpsro_data%sclf(nobs))
118 allocate(gpsro_data%record_number(nobs))
119 allocate(gpsro_data%bnd_obs(nobs))
120 allocate(gpsro_data%bnd_preqc(nobs))
121 allocate(gpsro_data%bnd_obserr(nobs))
122 allocate(gpsro_data%bnd_obserr_final(nobs))
123 allocate(gpsro_data%bnd_obserr_adjust(nobs))
124 allocate(gpsro_data%bnd_gsihofx(nobs))
125 
126 call check(nf90_inq_varid(ncid_in, 'record_number@MetaData', varid_meta))
127 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%record_number, (/1/),(/nobs/)))
128 call check(nf90_inq_varid(ncid_in, 'latitude@MetaData', varid_meta))
129 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%lat(:), (/1/),(/nobs/)))
130 call check(nf90_inq_varid(ncid_in, 'longitude@MetaData', varid_meta))
131 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%lon(:), (/1/),(/nobs/)))
132 call check(nf90_inq_varid(ncid_in, 'time@MetaData', varid_meta))
133 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%time(:), (/1/),(/nobs/)))
134 call check(nf90_inq_varid(ncid_in, 'impact_height@MetaData', varid_meta))
135 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%impact_height, (/1/),(/nobs/)))
136 call check(nf90_inq_varid(ncid_in, 'impact_parameter@MetaData', varid_meta))
137 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%impact_parameter, (/1/),(/nobs/)))
138 call check(nf90_inq_varid(ncid_in, 'earth_radius_of_curvature@MetaData',varid_meta))
139 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%rfict, (/1/), (/nobs/)))
140 call check(nf90_inq_varid(ncid_in,'geoid_height_above_reference_ellipsoid@MetaData', varid_meta))
141 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%geoid, (/1/), (/nobs/)))
142 call check(nf90_inq_varid(ncid_in, 'occulting_sat_id@MetaData', varid_meta))
143 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%said, (/1/), (/nobs/)))
144 call check(nf90_inq_varid(ncid_in, 'occulting_sat_is@MetaData', varid_meta))
145 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%siid, (/1/), (/nobs/)))
146 call check(nf90_inq_varid(ncid_in, 'process_center@MetaData', varid_meta))
147 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%ogce, (/1/), (/nobs/)))
148 call check(nf90_inq_varid(ncid_in, 'reference_sat_id', varid_meta))
149 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%ptid, (/1/), (/nobs/)))
150 call check(nf90_inq_varid(ncid_in, 'gnss_sat_class@MetaData', varid_meta))
151 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%sclf, (/1/), (/nobs/)))
152 call check(nf90_inq_varid(ncid_in, 'sensor_azimuth_angle@MetaData',varid_meta))
153 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%azim, (/1/), (/nobs/)))
154 call check(nf90_inq_varid(ncid_in, 'ascending_flag@MetaData', varid_meta))
155 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%asce, (/1/), (/nobs/)))
156 call check(nf90_inq_varid(ncid_in, 'bending_angle@ObsValue', varid_meta))
157 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%bnd_obs, (/1/),(/nobs/)))
158 call check(nf90_inq_varid(ncid_in, 'bending_angle@PreQC', varid_meta))
159 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%bnd_preqc, (/1/),(/nobs/)))
160 call check(nf90_inq_varid(ncid_in, 'bending_angle@GsiFinalObsErrorInv',varid_meta))
161 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%bnd_obserr_final, (/1/),(/nobs/)))
162 call check(nf90_inq_varid(ncid_in, 'bending_angle@GsiAdjustObsErrorInv',varid_meta))
163 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%bnd_obserr_adjust, (/1/),(/nobs/)))
164 call check(nf90_inq_varid(ncid_in, 'bending_angle@GsiHofX', varid_meta))
165 call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%bnd_gsihofx, (/1/),(/nobs/)))
166 
167 where(gpsro_data%bnd_obserr_adjust>0.0 )
168  gpsro_data%bnd_obserr_adjust = 1.0/gpsro_data%bnd_obserr_adjust
169 elsewhere
170  gpsro_data%bnd_obserr_adjust = huge_single
171 endwhere
172 
173 where(gpsro_data%bnd_obserr_final>0.0 )
174  gpsro_data%bnd_obserr_final = 1.0/gpsro_data%bnd_obserr_final
175 elsewhere
176  gpsro_data%bnd_obserr_final = huge_single
177 endwhere
178 gpsro_data%bnd_obserr = gpsro_data%bnd_obserr_adjust
179 
180 if ( geovalwrite == "1" ) then
181  geoname_out='./gnssro_geoval_'//trim(anatime)//'.nc4'
182  allocate(gpsro_data%air_temperature(nlevs,nobs))
183  allocate(gpsro_data%virtual_temperature(nlevs,nobs))
184  allocate(gpsro_data%specific_humidity(nlevs,nobs))
185  allocate(gpsro_data%geopotential_height(nlevs,nobs))
186  allocate(gpsro_data%geopotential_height_levels(nlevs1,nobs))
187  allocate(gpsro_data%air_pressure(nlevs,nobs))
188  allocate(gpsro_data%air_pressure_levels(nlevs1,nobs))
189  allocate(gpsro_data%surface_geopotential_height(nobs))
190  call check(nf90_inq_varid(ncid_in, 'surface_geopotential_height', varid_meta))
191  call check(nf90_get_var(ncid_in, varid_meta,gpsro_data%surface_geopotential_height, (/1/), (/nobs/)))
192  call check(nf90_inq_varid(ncid_in, 'air_temperature', varid_geoval_nlevs))
193  call check(nf90_get_var(ncid_in, varid_geoval_nlevs,gpsro_data%air_temperature, start=(/1,1/), count=(/nlevs,nobs/)))
194  call check(nf90_inq_varid(ncid_in, 'virtual_temperature', varid_geoval_nlevs))
195  call check(nf90_get_var(ncid_in,varid_geoval_nlevs,gpsro_data%virtual_temperature, start=(/1,1/),count=(/nlevs,nobs/)))
196  call check(nf90_inq_varid(ncid_in, 'specific_humidity', varid_meta))
197  call check(nf90_get_var(ncid_in, varid_meta, gpsro_data%specific_humidity,start=(/1,1/), count=(/nlevs,nobs/)))
198  call check(nf90_inq_varid(ncid_in, 'air_pressure', varid_geoval_nlevs))
199  call check(nf90_get_var(ncid_in, varid_geoval_nlevs, gpsro_data%air_pressure,start=(/1,1/), count=(/nlevs, nobs/)))
200  call check(nf90_inq_varid(ncid_in, 'air_pressure_levels', varid_geoval_nlevs1))
201  call check(nf90_get_var(ncid_in, varid_geoval_nlevs1,gpsro_data%air_pressure_levels, start=(/1,1/), count=(/nlevs1, nobs/)))
202  call check(nf90_inq_varid(ncid_in, 'geopotential_height', varid_geoval_nlevs))
203  call check(nf90_get_var(ncid_in, varid_geoval_nlevs,gpsro_data%geopotential_height, start=(/1,1/), count=(/nlevs,nobs/)))
204  call check(nf90_inq_varid(ncid_in, 'geopotential_height_levels',varid_geoval_nlevs1))
205  call check(nf90_get_var(ncid_in, varid_geoval_nlevs1,gpsro_data%geopotential_height_levels, start=(/1,1/), &
206  count=(/nlevs1, nobs/)))
207 
208  allocate(gpsro_subset%air_temperature(nlevs,nobs))
209  allocate(gpsro_subset%virtual_temperature(nlevs,nobs))
210  allocate(gpsro_subset%specific_humidity(nlevs,nobs))
211  allocate(gpsro_subset%geopotential_height(nlevs,nobs))
212  allocate(gpsro_subset%geopotential_height_levels(nlevs1,nobs))
213  allocate(gpsro_subset%air_pressure(nlevs,nobs))
214  allocate(gpsro_subset%air_pressure_levels(nlevs1,nobs))
215  allocate(gpsro_subset%surface_geopotential_height(nobs))
216 
217 end if
218 
219 allocate(gpsro_subset%lat(nobs))
220 allocate(gpsro_subset%lon(nobs))
221 allocate(gpsro_subset%time(nobs))
222 allocate(gpsro_subset%asce(nobs))
223 allocate(gpsro_subset%azim(nobs))
224 allocate(gpsro_subset%geoid(nobs))
225 allocate(gpsro_subset%rfict(nobs))
226 allocate(gpsro_subset%impact_height(nobs))
227 allocate(gpsro_subset%impact_parameter(nobs))
228 allocate(gpsro_subset%said(nobs))
229 allocate(gpsro_subset%siid(nobs))
230 allocate(gpsro_subset%ogce(nobs))
231 allocate(gpsro_subset%ptid(nobs))
232 allocate(gpsro_subset%sclf(nobs))
233 allocate(gpsro_subset%record_number(nobs))
234 allocate(gpsro_subset%bnd_obs(nobs))
235 allocate(gpsro_subset%bnd_preqc(nobs))
236 allocate(gpsro_subset%bnd_obserr(nobs))
237 allocate(gpsro_subset%bnd_obserr_adjust(nobs))
238 allocate(gpsro_subset%bnd_obserr_final(nobs))
239 allocate(gpsro_subset%bnd_gsihofx(nobs))
240 
241 nsub = 0
242 nrecs = maxval(gpsro_data%record_number(:))
243 
244 do j = 1, maxval(gpsro_data%record_number(:))
245  do i = 1, nobs
246  if ( gpsro_data%record_number(i) .eq. j ) then
247  nsub=nsub+1
248  gpsro_subset%lat(nsub) = gpsro_data%lat(i)
249  gpsro_subset%lon(nsub) = gpsro_data%lon(i)
250  gpsro_subset%time(nsub) = gpsro_data%time(i)
251  if (gpsro_subset%time(nsub) .eq. -3.00) gpsro_subset%time(nsub) = -2.999
252  gpsro_subset%impact_height(nsub) = gpsro_data%impact_height(i)
253  gpsro_subset%impact_parameter(nsub) = gpsro_data%impact_parameter(i)
254  gpsro_subset%geoid(nsub) = gpsro_data%geoid(i)
255  gpsro_subset%rfict(nsub) = gpsro_data%rfict(i)
256  gpsro_subset%said(nsub) = gpsro_data%said(i)
257  gpsro_subset%siid(nsub) = gpsro_data%siid(i)
258  gpsro_subset%ogce(nsub) = gpsro_data%ogce(i)
259  gpsro_subset%ptid(nsub) = gpsro_data%ptid(i)
260  gpsro_subset%sclf(nsub) = gpsro_data%sclf(i)
261  gpsro_subset%record_number(nsub) = gpsro_data%record_number(i)
262  gpsro_subset%azim(nsub) = gpsro_data%azim(i)
263  gpsro_subset%asce(nsub) = gpsro_data%asce(i)
264  gpsro_subset%bnd_obs(nsub) = gpsro_data%bnd_obs(i)
265  gpsro_subset%bnd_preqc(nsub) = gpsro_data%bnd_preqc(i)
266  gpsro_subset%bnd_obserr(nsub) = gpsro_data%bnd_obserr(i)
267  gpsro_subset%bnd_obserr_final(nsub) =gpsro_data%bnd_obserr_final(i)
268  gpsro_subset%bnd_obserr_adjust(nsub) = gpsro_data%bnd_obserr_adjust(i)
269  gpsro_subset%bnd_gsihofx(nsub) = gpsro_data%bnd_gsihofx(i)
270  if ( geovalwrite == "1" ) then
271  gpsro_subset%surface_geopotential_height(nsub) = gpsro_data%surface_geopotential_height(i)
272  gpsro_subset%air_temperature(:,nsub) = gpsro_data%air_temperature(:,i)
273  gpsro_subset%virtual_temperature(:,nsub) = gpsro_data%virtual_temperature(:,i)
274  gpsro_subset%specific_humidity(:,nsub) = gpsro_data%specific_humidity(:, i )
275  gpsro_subset%air_pressure(:,nsub) = gpsro_data%air_pressure(:, i )
276  gpsro_subset%air_pressure_levels(:,nsub) = gpsro_data%air_pressure_levels(:, i )
277  gpsro_subset%geopotential_height(:,nsub) = gpsro_data%geopotential_height(:, i )
278  gpsro_subset%geopotential_height_levels(:,nsub) = gpsro_data%geopotential_height_levels(:,i)
279  endif
280 
281  end if
282  end do
283 end do
284 
285 call check( nf90_create(trim(obsname_out), nf90_clobber, ncid_out))
286 call check( nf90_put_att(ncid_out, nf90_global, 'date_time', anatime_i) )
287 call check( nf90_def_dim(ncid_out, 'nlocs', nsub, nlocs_dimid) )
288 call check( nf90_def_var(ncid_out, "latitude@MetaData", nf90_float, nlocs_dimid, varid_lat) )
289 call check( nf90_def_var(ncid_out, "longitude@MetaData", nf90_float, nlocs_dimid, varid_lon) )
290 call check( nf90_def_var(ncid_out, "time@MetaData", nf90_float, nlocs_dimid, varid_time) )
291 call check( nf90_def_var(ncid_out, "impact_height@MetaData", nf90_float, nlocs_dimid, varid_imph))
292 call check( nf90_put_att(ncid_out, varid_imph, "longname", "distance from mean sea level" ))
293 call check( nf90_put_att(ncid_out, varid_imph, "units", "Meters" ))
294 call check( nf90_put_att(ncid_out, varid_imph, "valid_range", "0 - 200 km" ))
295 call check( nf90_def_var(ncid_out, "impact_parameter@MetaData",nf90_float, nlocs_dimid, varid_impp))
296 call check( nf90_put_att(ncid_out, varid_impp, "longname", "distance from centre of curvature" ))
297 call check( nf90_put_att(ncid_out, varid_impp, "units", "Meters" ))
298 call check( nf90_put_att(ncid_out, varid_impp, "valid_range", "6200 - 6600 km" ))
299 call check( nf90_def_var(ncid_out, "geoid_height_above_reference_ellipsoid@MetaData", nf90_float, nlocs_dimid, varid_geoid))
300 call check( nf90_put_att(ncid_out, varid_geoid, "longname", "Geoid height above WGS-84 ellipsoid" ))
301 call check( nf90_put_att(ncid_out, varid_geoid, "units", "Meters" ))
302 call check( nf90_put_att(ncid_out, varid_geoid, "valid_range", "-200 - 200 m" ))
303 call check( nf90_def_var(ncid_out, "earth_radius_of_curvature@MetaData", nf90_float, nlocs_dimid, varid_rfict))
304 call check( nf90_put_att(ncid_out, varid_rfict, "longname", ’"Earths local radius of curvature" ))
305 call check( nf90_put_att(ncid_out, varid_rfict, "units", "Meters" ))
306 call check( nf90_put_att(ncid_out, varid_rfict, "valid_range", "6200 - 6600 km" ))
307 call check( nf90_def_var(ncid_out, "sensor_azimuth_angle@MetaData",nf90_float, nlocs_dimid, varid_azim))
308 call check( nf90_put_att(ncid_out, varid_azim, "longname", "GNSS->LEO line of sight" ))
309 call check( nf90_put_att(ncid_out, varid_azim, "units", "Degree" ))
310 call check( nf90_put_att(ncid_out, varid_azim, "valid_range", "0 - 360 degree" ))
311 call check( nf90_def_var(ncid_out, "record_number@MetaData", nf90_int,nlocs_dimid, varid_rec))
312 call check( nf90_put_att(ncid_out, varid_rec, "longname", "GNSS RO profile identifier" ))
313 call check( nf90_def_var(ncid_out, "occulting_sat_id@MetaData", nf90_int, nlocs_dimid, varid_said))
314 call check( nf90_put_att(ncid_out, varid_said, "longname", "Low Earth Orbit satellite identifier, e.g., COSMIC2=750-755" ))
315 call check( nf90_def_var(ncid_out, "occulting_sat_is@MetaData", nf90_int, nlocs_dimid, varid_siid))
316 call check( nf90_put_att(ncid_out, varid_siid, "longname", "satellite instrument" ))
317 call check( nf90_def_var(ncid_out, "process_center@MetaData", nf90_int, nlocs_dimid, varid_ogce))
318 call check( nf90_put_att(ncid_out, varid_ogce, "longname", "originally data processing_center, e.g., 60 for UCAR, 94 for DMI" ))
319 call check( nf90_def_var(ncid_out, "reference_sat_id@MetaData", nf90_int, nlocs_dimid, varid_ptid))
320 call check( nf90_put_att(ncid_out, varid_ptid, "longname", "GNSS satellite transmitter identifier (1-32)" ))
321 call check( nf90_def_var(ncid_out, "gnss_sat_class@MetaData", nf90_int, nlocs_dimid, varid_sclf))
322 call check( nf90_put_att(ncid_out, varid_sclf, "longname", "GNSS satellite classification, e.g, 401=GPS, 402=GLONASS" ))
323 call check( nf90_def_var(ncid_out, "ascending_flag@MetaData", nf90_int, nlocs_dimid, varid_asce))
324 call check( nf90_put_att(ncid_out, varid_asce, "longname", "the original occultation ascending/descending flag" ))
325 call check( nf90_def_var(ncid_out, "bending_angle@ObsValue", nf90_float, nlocs_dimid, varid_bnd))
326 call check( nf90_put_att(ncid_out, varid_bnd, "longname", "Bending Angle" ))
327 call check( nf90_put_att(ncid_out, varid_bnd, "units", "Radians" ))
328 call check( nf90_put_att(ncid_out, varid_bnd, "valid_range", "-0.001 - 0.08 Radians" ))
329 call check( nf90_def_var(ncid_out, "bending_angle@ObsError", nf90_float, nlocs_dimid, varid_obserr))
330 call check( nf90_put_att(ncid_out, varid_obserr, "longname", "GSI final observation error" ))
331 call check( nf90_put_att(ncid_out, varid_obserr, "units", "Radians" ))
332 call check( nf90_put_att(ncid_out, varid_obserr, "valid_range", "0 - 0.008 Radians" ))
333 call check( nf90_def_var(ncid_out, "bending_angle@GsiAdjustObsError", nf90_float, nlocs_dimid, varid_obserr_adjust))
334 call check( nf90_put_att(ncid_out, varid_obserr_adjust, "longname", "GSI adjust observation error" ))
335 call check( nf90_put_att(ncid_out, varid_obserr_adjust, "units", "Radians" ))
336 call check( nf90_put_att(ncid_out, varid_obserr_adjust, "valid_range", "0 - 0.008 Radians" ))
337 call check( nf90_def_var(ncid_out, "bending_angle@GsiFinalObsError", nf90_float, nlocs_dimid, varid_obserr_final))
338 call check( nf90_put_att(ncid_out, varid_obserr_final, "longname", "GSI final observation error (inflated)" ))
339 call check( nf90_put_att(ncid_out, varid_obserr_final, "units", "Radians" ))
340 call check( nf90_put_att(ncid_out, varid_obserr_final, "valid_range", "0 - 0.008 Radians" ))
341 call check( nf90_def_var(ncid_out, "bending_angle@PreQC", nf90_int,nlocs_dimid, varid_preqc))
342 call check( nf90_put_att(ncid_out, varid_preqc, "longname", "GSI output QC flags" ))
343 call check( nf90_put_att(ncid_out, varid_preqc, "valid_range", "1 - 5" ))
344 call check( nf90_def_var(ncid_out, "bending_angle@GsiHofX", nf90_float, nlocs_dimid, varid_gsihofx))
345 call check( nf90_enddef(ncid_out) )
346 
347 if ( geovalwrite =="1" ) then
348 call check( nf90_create(trim(geoname_out), nf90_clobber, ncid_out2))
349 call check( nf90_put_att(ncid_out2, nf90_global, 'date_time', anatime_i) )
350 call check( nf90_def_dim(ncid_out2, 'nlocs', nsub, nlocs_dimid2) )
351 call check( nf90_def_dim(ncid_out2, 'nlevs', nlevs, nlevs_dimid) )
352 call check( nf90_def_dim(ncid_out2, 'nlevs1', nlevs1, nlevs1_dimid) )
353 call check( nf90_def_var(ncid_out2, "latitude", nf90_float, nlocs_dimid, varid_lat2) )
354 call check( nf90_def_var(ncid_out2, "longitude", nf90_float, nlocs_dimid, varid_lon2) )
355 call check( nf90_def_var(ncid_out2, "time", nf90_float, nlocs_dimid, varid_time2) )
356 call check( nf90_def_var(ncid_out2, "surface_geopotential_height", nf90_float, nlocs_dimid2, varid_geo_geop_sfc) )
357 call check( nf90_def_var(ncid_out2, "air_temperature", nf90_float, (/nlevs_dimid,nlocs_dimid2/), varid_geo_airtemp) )
358 call check( nf90_def_var(ncid_out2, "virtual_temperature", nf90_float, (/nlevs_dimid,nlocs_dimid2/), varid_geo_virtemp) )
359 call check( nf90_def_var(ncid_out2, "specific_humidity", nf90_float, (/nlevs_dimid,nlocs_dimid2/), varid_geo_shum) )
360 call check( nf90_def_var(ncid_out2, "air_pressure", nf90_float, (/nlevs_dimid,nlocs_dimid2/), varid_geo_pres))
361 call check( nf90_def_var(ncid_out2, "air_pressure_levels", nf90_float, (/nlevs1_dimid,nlocs_dimid2/), varid_geo_pres1))
362 call check( nf90_def_var(ncid_out2, "geopotential_height", nf90_float, (/nlevs_dimid,nlocs_dimid2/), varid_geo_geop))
363 call check( nf90_def_var(ncid_out2, "geopotential_height_levels", nf90_float, (/nlevs1_dimid,nlocs_dimid2/), varid_geo_geop1))
364 call check( nf90_enddef(ncid_out2) )
365 end if
366 
367 call check( nf90_put_var(ncid_out, varid_lat, gpsro_subset%lat(1:nsub)) )
368 call check( nf90_put_var(ncid_out, varid_lon, gpsro_subset%lon(1:nsub)) )
369 call check( nf90_put_var(ncid_out, varid_time, gpsro_subset%time(1:nsub)) )
370 call check( nf90_put_var(ncid_out, varid_imph, gpsro_subset%impact_height(1:nsub)) )
371 call check( nf90_put_var(ncid_out, varid_impp, gpsro_subset%impact_parameter(1:nsub)) )
372 call check( nf90_put_var(ncid_out, varid_geoid, gpsro_subset%geoid(1:nsub)) )
373 call check( nf90_put_var(ncid_out, varid_rfict, gpsro_subset%rfict(1:nsub)) )
374 call check( nf90_put_var(ncid_out, varid_azim, gpsro_subset%azim(1:nsub)) )
375 call check( nf90_put_var(ncid_out, varid_rec, gpsro_subset%record_number(1:nsub)) )
376 call check( nf90_put_var(ncid_out, varid_said, gpsro_subset%said(1:nsub)) )
377 call check( nf90_put_var(ncid_out, varid_siid, gpsro_subset%siid(1:nsub)) )
378 call check( nf90_put_var(ncid_out, varid_ogce, gpsro_subset%ogce(1:nsub)) )
379 call check( nf90_put_var(ncid_out, varid_ptid, gpsro_subset%ptid(1:nsub)) )
380 call check( nf90_put_var(ncid_out, varid_sclf, gpsro_subset%sclf(1:nsub)) )
381 call check( nf90_put_var(ncid_out, varid_asce, gpsro_subset%asce(1:nsub)) )
382 call check( nf90_put_var(ncid_out, varid_bnd, gpsro_subset%bnd_obs(1:nsub)) )
383 call check( nf90_put_var(ncid_out, varid_preqc, gpsro_subset%bnd_preqc(1:nsub)) )
384 call check( nf90_put_var(ncid_out, varid_obserr, gpsro_subset%bnd_obserr(1:nsub)) )
385 call check( nf90_put_var(ncid_out, varid_obserr_final, gpsro_subset%bnd_obserr_final(1:nsub)) )
386 call check( nf90_put_var(ncid_out, varid_obserr_adjust, gpsro_subset%bnd_obserr_adjust(1:nsub)) )
387 call check( nf90_put_var(ncid_out, varid_gsihofx,gpsro_subset%bnd_gsihofx(1:nsub)) )
388 call check( nf90_close(ncid_out) )
389 
390 if ( geovalwrite == "1" ) then
391  call check( nf90_put_var(ncid_out2, varid_lat2, gpsro_subset%lat(1:nsub)) )
392  call check( nf90_put_var(ncid_out2, varid_lon2, gpsro_subset%lon(1:nsub)) )
393  call check( nf90_put_var(ncid_out2, varid_time2, gpsro_subset%time(1:nsub)) )
394  call check( nf90_put_var(ncid_out2, varid_geo_geop_sfc,gpsro_subset%surface_geopotential_height(1:nsub)) )
395  call check( nf90_put_var(ncid_out2, varid_geo_airtemp,gpsro_subset%air_temperature(1:nlevs,1:nsub)) )
396  call check( nf90_put_var(ncid_out2, varid_geo_virtemp,gpsro_subset%virtual_temperature(1:nlevs,1:nsub)) )
397  call check( nf90_put_var(ncid_out2, varid_geo_shum,gpsro_subset%specific_humidity(1:nlevs,1:nsub)) )
398  call check( nf90_put_var(ncid_out2, varid_geo_pres,gpsro_subset%air_pressure(1:nlevs,1:nsub)) )
399  call check( nf90_put_var(ncid_out2, varid_geo_geop,gpsro_subset%geopotential_height(1:nlevs,1:nsub)) )
400  call check( nf90_put_var(ncid_out2, varid_geo_pres1,gpsro_subset%air_pressure_levels(1:nlevs1,1:nsub)) )
401  call check( nf90_put_var(ncid_out2, varid_geo_geop1,gpsro_subset%geopotential_height_levels(1:nlevs1,1:nsub)) )
402  call check( nf90_close(ncid_out2) )
403 end if
404 
405 deallocate(gpsro_data%lat)
406 deallocate(gpsro_data%lon)
407 deallocate(gpsro_data%time)
408 deallocate(gpsro_data%asce)
409 deallocate(gpsro_data%azim)
410 deallocate(gpsro_data%geoid)
411 deallocate(gpsro_data%rfict)
412 deallocate(gpsro_data%impact_height)
413 deallocate(gpsro_data%impact_parameter)
414 deallocate(gpsro_data%said)
415 deallocate(gpsro_data%siid)
416 deallocate(gpsro_data%ogce)
417 deallocate(gpsro_data%ptid)
418 deallocate(gpsro_data%sclf)
419 deallocate(gpsro_data%record_number)
420 deallocate(gpsro_data%bnd_obs)
421 deallocate(gpsro_data%bnd_preqc)
422 deallocate(gpsro_data%bnd_obserr)
423 deallocate(gpsro_data%bnd_obserr_final)
424 deallocate(gpsro_data%bnd_obserr_adjust)
425 deallocate(gpsro_data%bnd_gsihofx)
426 
427 deallocate(gpsro_subset%lat)
428 deallocate(gpsro_subset%lon)
429 deallocate(gpsro_subset%time)
430 deallocate(gpsro_subset%asce)
431 deallocate(gpsro_subset%azim)
432 deallocate(gpsro_subset%geoid)
433 deallocate(gpsro_subset%rfict)
434 deallocate(gpsro_subset%impact_height)
435 deallocate(gpsro_subset%impact_parameter)
436 deallocate(gpsro_subset%said)
437 deallocate(gpsro_subset%siid)
438 deallocate(gpsro_subset%ogce)
439 deallocate(gpsro_subset%ptid)
440 deallocate(gpsro_subset%sclf)
441 deallocate(gpsro_subset%record_number)
442 deallocate(gpsro_subset%bnd_obs)
443 deallocate(gpsro_subset%bnd_preqc)
444 deallocate(gpsro_subset%bnd_obserr)
445 deallocate(gpsro_subset%bnd_obserr_adjust)
446 deallocate(gpsro_subset%bnd_obserr_final)
447 deallocate(gpsro_subset%bnd_gsihofx)
448 
449 if ( geovalwrite == "1" ) then
450  deallocate(gpsro_data%air_temperature)
451  deallocate(gpsro_data%virtual_temperature)
452  deallocate(gpsro_data%specific_humidity)
453  deallocate(gpsro_data%geopotential_height)
454  deallocate(gpsro_data%geopotential_height_levels)
455  deallocate(gpsro_data%air_pressure)
456  deallocate(gpsro_data%air_pressure_levels)
457  deallocate(gpsro_data%surface_geopotential_height)
458  deallocate(gpsro_subset%air_temperature)
459  deallocate(gpsro_subset%virtual_temperature)
460  deallocate(gpsro_subset%specific_humidity)
461  deallocate(gpsro_subset%geopotential_height)
462  deallocate(gpsro_subset%geopotential_height_levels)
463  deallocate(gpsro_subset%air_pressure)
464  deallocate(gpsro_subset%air_pressure_levels)
465  deallocate(gpsro_subset%surface_geopotential_height)
466 end if
467 
468 contains
469  subroutine check(status)
470  integer, intent ( in) :: status
471 
472  if(status /= nf90_noerr) then
473  print *, trim(nf90_strerror(status))
474  stop "Stopped"
475  end if
476 
477  end subroutine check
478 
479 
480 end program gnssro_gsidiag_full2small
count
static void count(void *counter, const double *data, size_t n)
Definition: UnitTests.cc:531
check
subroutine check(status)
Definition: gnssro_bufr2ioda.f90:443
gnssro_gsidiag_full2small
program gnssro_gsidiag_full2small
Definition: gnssro_gsidiag2ioda.f90:17
gpsro_type
Definition: gnssro_bufr2ioda.f90:57