FV3DYCORE  Version1.0.0
sim_nc_mod.F90
Go to the documentation of this file.
1 !***********************************************************************
2 !* GNU Lesser General Public License
3 !*
4 !* This file is part of the FV3 dynamical core.
5 !*
6 !* The FV3 dynamical core is free software: you can redistribute it
7 !* and/or modify it under the terms of the
8 !* GNU Lesser General Public License as published by the
9 !* Free Software Foundation, either version 3 of the License, or
10 !* (at your option) any later version.
11 !*
12 !* The FV3 dynamical core is distributed in the hope that it will be
13 !* useful, but WITHOUT ANYWARRANTY; without even the implied warranty
14 !* of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
15 !* See the GNU General Public License for more details.
16 !*
17 !* You should have received a copy of the GNU Lesser General Public
18 !* License along with the FV3 dynamical core.
19 !* If not, see <http://www.gnu.org/licenses/>.
20 !***********************************************************************
21 
27 
28 module sim_nc_mod
29 
30 ! This is S-J Lin's private netcdf file reader
31 ! This code is needed because FMS utility (read_data) led to too much
32 ! memory usage and too many files openned. Perhaps lower-level FMS IO
33 ! calls should be used instead.
34 
35 #if defined(OLD_PT_TO_T) || defined(OLD_COS_SG)
36 #error
37 #error Compile time options -DOLD_PT_TO_T and -DOLD_COS_SG are no longer supported. Please remove them from your XML.
38 #error
39 #endif
40 
41  use mpp_mod, only: mpp_error, fatal
42 
43  implicit none
44 #include <netcdf.inc>
45 
46  private
51 
52  contains
53 
54  subroutine open_ncfile( iflnm, ncid )
55  character(len=*), intent(in):: iflnm
56  integer, intent(out):: ncid
57  integer:: status
58 
59  status = nf_open(iflnm, nf_nowrite, ncid)
60  if (status .ne. nf_noerr) call handle_err('nf_open',status)
61 
62 
63  end subroutine open_ncfile
64 
65 
66  subroutine close_ncfile( ncid )
67  integer, intent(in):: ncid
68  integer:: status
69 
70  status = nf_close(ncid)
71  if (status .ne. nf_noerr) call handle_err('nf_close',status)
72 
73 
74  end subroutine close_ncfile
75 
76 
77  subroutine get_ncdim1( ncid, var1_name, im )
78  integer, intent(in):: ncid
79  character(len=*), intent(in):: var1_name
80  integer, intent(out):: im
81  integer:: status, var1id
82 
83  status = nf_inq_dimid(ncid, var1_name, var1id)
84  if (status .ne. nf_noerr) call handle_err('dimid '//var1_name,status)
85 
86  status = nf_inq_dimlen(ncid, var1id, im)
87  if (status .ne. nf_noerr) call handle_err('dimid '//var1_name,status)
88 
89  end subroutine get_ncdim1
90 
92  subroutine get_var1_double( ncid, var1_name, im, var1, var_exist )
93  integer, intent(in):: ncid
94  character(len=*), intent(in):: var1_name
95  integer, intent(in):: im
96  logical, intent(out), optional:: var_exist
97  real(kind=8), intent(out):: var1(im)
98  integer:: status, var1id
99 
100  status = nf_inq_varid(ncid, var1_name, var1id)
101  if (status .ne. nf_noerr) then
102 ! call handle_err('varid '//var1_name,status)
103  if(present(var_exist) ) var_exist = .false.
104  else
105  status = nf_get_var_double(ncid, var1id, var1)
106  if (status .ne. nf_noerr) call handle_err('varid '//var1_name,status)
107  if(present(var_exist) ) var_exist = .true.
108  endif
109 
110 
111  end subroutine get_var1_double
112 
113 
114 ! 4-byte data:
115  subroutine get_var1_real( ncid, var1_name, im, var1, var_exist )
116  integer, intent(in):: ncid
117  character(len=*), intent(in):: var1_name
118  integer, intent(in):: im
119  logical, intent(out), optional:: var_exist
120  real(kind=4), intent(out):: var1(im)
121  integer:: status, var1id
122 
123  status = nf_inq_varid(ncid, var1_name, var1id)
124  if (status .ne. nf_noerr) then
125 ! call handle_err(status)
126  if(present(var_exist) ) var_exist = .false.
127  else
128  status = nf_get_var_real(ncid, var1id, var1)
129  if (status .ne. nf_noerr) call handle_err('get_var1_real1 '//var1_name,status)
130  if(present(var_exist) ) var_exist = .true.
131  endif
132 
133 
134  end subroutine get_var1_real
135 
136  subroutine get_var2_real( ncid, var_name, im, jm, var2 )
137  integer, intent(in):: ncid
138  character(len=*), intent(in):: var_name
139  integer, intent(in):: im, jm
140  real(kind=4), intent(out):: var2(im)
141 
142  integer:: status, var1id
143 
144  status = nf_inq_varid(ncid, var_name, var1id)
145  if (status .ne. nf_noerr) call handle_err('get_var2_real varid '//var_name,status)
146 
147  status = nf_get_var_real(ncid, var1id, var2)
148  if (status .ne. nf_noerr) call handle_err('get_var2_real get_var'//var_name,status)
149 
150  end subroutine get_var2_real
151 
152  subroutine get_var2_r4( ncid, var2_name, is,ie, js,je, var2, time_slice )
153  integer, intent(in):: ncid
154  character(len=*), intent(in):: var2_name
155  integer, intent(in):: is, ie, js, je
156  real(kind=4), intent(out):: var2(is:ie,js:je)
157  integer, intent(in), optional :: time_slice
158 !
159  real(kind=4), dimension(1) :: time
160  integer, dimension(3):: start, nreco
161  integer:: status, var2id
162 
163  status = nf_inq_varid(ncid, var2_name, var2id)
164  if (status .ne. nf_noerr) call handle_err('get_var2_r4 varid'//var2_name,status)
165 
166  start(1) = is; start(2) = js; start(3) = 1
167  if ( present(time_slice) ) then
168  start(3) = time_slice
169  end if
170 
171  nreco(1) = ie - is + 1
172  nreco(2) = je - js + 1
173  nreco(3) = 1
174 
175  status = nf_get_vara_real(ncid, var2id, start, nreco, var2)
176  if (status .ne. nf_noerr) call handle_err('get_var2_r4 get_vara_real'//var2_name,status)
177 
178  end subroutine get_var2_r4
179 
180  subroutine get_var2_double( ncid, var2_name, im, jm, var2 )
181  integer, intent(in):: ncid
182  character(len=*), intent(in):: var2_name
183  integer, intent(in):: im, jm
184  real(kind=8), intent(out):: var2(im,jm)
185 
186  integer:: status, var2id
187 
188  status = nf_inq_varid(ncid, var2_name, var2id)
189  if (status .ne. nf_noerr) call handle_err('get_var2_double varid'//var2_name,status)
190 
191  status = nf_get_var_double(ncid, var2id, var2)
192  if (status .ne. nf_noerr) call handle_err('get_var2_double get_var_double'//var2_name,status)
193 
194 
195  end subroutine get_var2_double
196 
197 
198  subroutine get_var3_double( ncid, var3_name, im, jm, km, var3 )
199  integer, intent(in):: ncid
200  character(len=*), intent(in):: var3_name
201  integer, intent(in):: im, jm, km
202  real(kind=8), intent(out):: var3(im,jm,km)
203 
204  integer:: status, var3id
205 
206  status = nf_inq_varid(ncid, var3_name, var3id)
207 
208  if (status .ne. nf_noerr) &
209  call handle_err('get_var3_double varid '//var3_name,status)
210 
211  status = nf_get_var_double(ncid, var3id, var3)
212  if (status .ne. nf_noerr) &
213  call handle_err('get_var3_double get_vara_double '//var3_name,status)
214 
215  end subroutine get_var3_double
216 
217  subroutine get_var3_real( ncid, var3_name, im, jm, km, var3 )
218  integer, intent(in):: ncid
219  character(len=*), intent(in):: var3_name
220  integer, intent(in):: im, jm, km
221  real(kind=4), intent(out):: var3(im,jm,km)
222 
223  integer:: status, var3id
224 
225  status = nf_inq_varid(ncid, var3_name, var3id)
226 
227  if (status .ne. nf_noerr) &
228  call handle_err('get_var3_real varid '//var3_name,status)
229  status = nf_get_var_real(ncid, var3id, var3)
230 
231  if (status .ne. nf_noerr) &
232  call handle_err('get_var3_real get_var_real '//var3_name,status)
233 
234  end subroutine get_var3_real
235 
236 
237  subroutine check_var_exists(ncid, var_name, status)
238  integer, intent(in):: ncid
239  integer, intent(inout) :: status
240  character(len=*), intent(in):: var_name
241  integer:: varid
242  status = nf_inq_varid(ncid, var_name, varid)
243  end subroutine check_var_exists
244 
245  subroutine get_var3_r4( ncid, var3_name, is,ie, js,je, ks,ke, var3, time_slice )
246  integer, intent(in):: ncid
247  character(len=*), intent(in):: var3_name
248  integer, intent(in):: is, ie, js, je, ks,ke
249  real(kind=4), intent(out):: var3(is:ie,js:je,ks:ke)
250  integer, intent(in), optional :: time_slice
251 !
252  real(kind=4), dimension(1) :: time
253  integer, dimension(4):: start, nreco
254  integer:: status, var3id
255 
256  status = nf_inq_varid(ncid, var3_name, var3id)
257  if (status .ne. nf_noerr) call handle_err('get_var3_r4 varid '//var3_name,status)
258 
259  start(1) = is; start(2) = js; start(3) = ks; start(4) = 1
260  if ( present(time_slice) ) then
261  start(4) = time_slice
262  end if
263 
264  nreco(1) = ie - is + 1
265  nreco(2) = je - js + 1
266  nreco(3) = ke - ks + 1
267  nreco(4) = 1
268 
269  status = nf_get_vara_real(ncid, var3id, start, nreco, var3)
270  if (status .ne. nf_noerr) call handle_err('get_var3_r4 get_vara_real '//var3_name,status)
271 
272  end subroutine get_var3_r4
273 
274 
275  subroutine get_var4_real( ncid, var4_name, im, jm, km, nt, var4 )
276  implicit none
277 #include <netcdf.inc>
278  integer, intent(in):: ncid
279  character*(*), intent(in):: var4_name
280  integer, intent(in):: im, jm, km, nt
281  real*4:: wk4(im,jm,km,4)
282  real*4, intent(out):: var4(im,jm)
283  integer:: status, var4id
284  integer:: start(4), icount(4)
285  integer:: i,j
286 
287  start(1) = 1
288  start(2) = 1
289  start(3) = 1
290  start(4) = nt
291 
292  icount(1) = im ! all range
293  icount(2) = jm ! all range
294  icount(3) = km ! all range
295  icount(4) = 1 ! one time level at a time
296 
297 ! write(*,*) nt, 'Within get_var4_double: ', var4_name
298 
299  status = nf_inq_varid(ncid, var4_name, var4id)
300 ! write(*,*) '#1', status, ncid, var4id
301 
302  status = nf_get_vara_real(ncid, var4id, start, icount, var4)
303 ! status = nf_get_vara_real(ncid, var4id, start, icount, wk4)
304 ! write(*,*) '#2', status, ncid, var4id
305 
306  do j=1,jm
307  do i=1,im
308 ! var4(i,j) = wk4(i,j,1,nt)
309  enddo
310  enddo
311 
312  if (status .ne. nf_noerr) call handle_err('get_var4_r4 get_vara_real '//var4_name,status)
313 
314  end subroutine get_var4_real
315 
316 
317  subroutine get_var4_double( ncid, var4_name, im, jm, km, nt, var4 )
318  integer, intent(in):: ncid
319  character(len=*), intent(in):: var4_name
320  integer, intent(in):: im, jm, km, nt
321  real(kind=8), intent(out):: var4(im,jm,km,1)
322  integer:: status, var4id
323 !
324  integer:: start(4), icount(4)
325 
326  start(1) = 1
327  start(2) = 1
328  start(3) = 1
329  start(4) = nt
330 
331  icount(1) = im ! all range
332  icount(2) = jm ! all range
333  icount(3) = km ! all range
334  icount(4) = 1 ! one time level at a time
335 
336  status = nf_inq_varid(ncid, var4_name, var4id)
337  status = nf_get_vara_double(ncid, var4id, start, icount, var4)
338 
339  if (status .ne. nf_noerr) call handle_err('get_var4_double get_vara_double '//var4_name,status)
340 
341  end subroutine get_var4_double
342 !------------------------------------------------------------------------
343 
344  subroutine get_real3( ncid, var4_name, im, jm, nt, var4 )
345 ! This is for multi-time-level 2D var
346  integer, intent(in):: ncid
347  character(len=*), intent(in):: var4_name
348  integer, intent(in):: im, jm, nt
349  real(kind=4), intent(out):: var4(im,jm)
350  integer:: status, var4id
351  integer:: start(3), icount(3)
352  integer:: i,j
353 
354  start(1) = 1
355  start(2) = 1
356  start(3) = nt
357 
358  icount(1) = im
359  icount(2) = jm
360  icount(3) = 1
361 
362  status = nf_inq_varid(ncid, var4_name, var4id)
363  status = nf_get_vara_real(ncid, var4id, start, icount, var4)
364 
365  if (status .ne. nf_noerr) &
366  call handle_err('get_real3 get_vara_real '//var4_name,status)
367 
368  end subroutine get_real3
369 !------------------------------------------------------------------------
370 
371  logical function check_var( ncid, var3_name)
372  integer, intent(in):: ncid
373  character(len=*), intent(in):: var3_name
374 
375  integer:: status, var3id
376 
377  status = nf_inq_varid(ncid, var3_name, var3id)
378  check_var = (status == nf_noerr)
379 
380  end function check_var
381 
382  subroutine get_var_att_str(ncid, var_name, att_name, att)
383  implicit none
384 #include <netcdf.inc>
385  integer, intent(in):: ncid
386  character*(*), intent(in):: var_name, att_name
387  character*(*), intent(out):: att
388 
389  integer:: status, varid
390 
391  status = nf_inq_varid(ncid, var_name, varid)
392  status = nf_get_att_text(ncid, varid, att_name, att)
393 
394  if (status .ne. nf_noerr) call handle_err('get_var_att_str '//var_name,status)
395 
396  end subroutine get_var_att_str
397 
398  subroutine get_var_att_double(ncid, var_name, att_name, value)
399  implicit none
400 #include <netcdf.inc>
401  integer, intent(in):: ncid
402  character*(*), intent(in):: var_name, att_name
403  real(kind=8), intent(out):: value
404 
405  integer:: status, varid
406 
407  status = nf_inq_varid(ncid, var_name, varid)
408  status = nf_get_att(ncid, varid, att_name, value)
409 
410  if (status .ne. nf_noerr) call handle_err('get_var_att_double '//var_name,status)
411 
412  end subroutine get_var_att_double
413 
414 
415  subroutine handle_err(idstr,status)
416  integer status
417  character(len=500) :: errstr
418  character(len=*) :: idstr
419 
420  if (status .ne. nf_noerr) then
421  write(errstr,*) 'Error in handle_err: ',trim(idstr)//' ',nf_strerror(status)
422  call mpp_error(fatal,errstr)
423  endif
424 
425  end subroutine handle_err
426 
428  subroutine calendar(year, month, day, hour)
429  integer, intent(inout) :: year ! year
430  integer, intent(inout) :: month ! month
431  integer, intent(inout) :: day ! day
432  integer, intent(inout) :: hour
433 !
434 ! Local variables
435 !
436  integer irem4,irem100
437  integer mdays(12)
438  data mdays /31,28,31,30,31,30,31,31,30,31,30,31/
439 !**** consider leap year
440 !
441  irem4 = mod( year, 4 )
442  irem100 = mod( year, 100 )
443  if( irem4 == 0 .and. irem100 /= 0) mdays(2) = 29
444 !
445  if( hour >= 24 ) then
446  day = day + 1
447  hour = hour - 24
448  end if
449 
450  if( day > mdays(month) ) then
451  day = day - mdays(month)
452  month = month + 1
453  end if
454  if( month > 12 ) then
455  year = year + 1
456  month = 1
457  end if
458 
459  end subroutine calendar
460 
461 end module sim_nc_mod
subroutine, public get_var3_r4(ncid, var3_name, is, ie, js, je, ks, ke, var3, time_slice)
Definition: sim_nc_mod.F90:246
subroutine, public handle_err(idstr, status)
Definition: sim_nc_mod.F90:416
subroutine, public get_var1_real(ncid, var1_name, im, var1, var_exist)
Definition: sim_nc_mod.F90:116
logical function, public check_var(ncid, var3_name)
Definition: sim_nc_mod.F90:372
subroutine, public get_var1_double(ncid, var1_name, im, var1, var_exist)
The &#39;get_var&#39; subroutines read in variables from netcdf files.
Definition: sim_nc_mod.F90:93
subroutine, public get_var2_double(ncid, var2_name, im, jm, var2)
Definition: sim_nc_mod.F90:181
subroutine get_var_att_str(ncid, var_name, att_name, att)
Definition: sim_nc_mod.F90:383
subroutine, public get_var3_real(ncid, var3_name, im, jm, km, var3)
Definition: sim_nc_mod.F90:218
The module &#39;sim_nc&#39; is a netcdf file reader.
Definition: sim_nc_mod.F90:28
subroutine, public get_var2_real(ncid, var_name, im, jm, var2)
Definition: sim_nc_mod.F90:137
subroutine get_var4_double(ncid, var4_name, im, jm, km, nt, var4)
Definition: sim_nc_mod.F90:318
subroutine, public get_var2_r4(ncid, var2_name, is, ie, js, je, var2, time_slice)
Definition: sim_nc_mod.F90:153
subroutine, public open_ncfile(iflnm, ncid)
Definition: sim_nc_mod.F90:55
subroutine, public get_var_att_double(ncid, var_name, att_name, value)
Definition: sim_nc_mod.F90:399
subroutine get_real3(ncid, var4_name, im, jm, nt, var4)
Definition: sim_nc_mod.F90:345
subroutine, public get_var3_double(ncid, var3_name, im, jm, km, var3)
Definition: sim_nc_mod.F90:199
subroutine, public get_ncdim1(ncid, var1_name, im)
Definition: sim_nc_mod.F90:78
subroutine, public check_var_exists(ncid, var_name, status)
Definition: sim_nc_mod.F90:238
subroutine calendar(year, month, day, hour)
The subroutine &#39;calendar&#39; computes the current GMT.
Definition: sim_nc_mod.F90:429
subroutine get_var4_real(ncid, var4_name, im, jm, km, nt, var4)
Definition: sim_nc_mod.F90:276
subroutine, public close_ncfile(ncid)
Definition: sim_nc_mod.F90:67