source: trunk/LMDZ.MARS/libf/phymars/iostart.F90 @ 3325

Last change on this file since 3325 was 2952, checked in by romain.vande, 19 months ago

Mars PCM + newstart + start_archive: Correct writing of variables inertiesoil and fluxgeo following r2919 and r2942
RV

File size: 37.1 KB
Line 
1MODULE iostart
2
3    IMPLICIT NONE
4PRIVATE
5    INTEGER,SAVE :: nid_start ! NetCDF file identifier for startfi.nc file
6    INTEGER,SAVE :: nid_restart ! NetCDF file identifier for restartfi.nc file
7    LOGICAL,SAVE :: first_time1d=.false.
8   
9    ! restartfi.nc file dimension identifiers: (see open_restartphy())
10    INTEGER,SAVE :: idim1 ! "index" dimension
11    INTEGER,SAVE :: idim2 ! "physical_points" dimension
12    INTEGER,SAVE :: idim3 ! "subsurface_layers" dimension
13    INTEGER,SAVE :: idim4 ! "nlayer_plus_1" dimension
14    INTEGER,SAVE :: idim5 ! "number_of_advected_fields" dimension
15    INTEGER,SAVE :: idim6 ! "nlayer" dimension
16    INTEGER,SAVE :: idim7 ! "Time" dimension
17    INTEGER,SAVE :: idim8 ! "nslope" dimension
18    INTEGER,SAVE :: idim9 ! "nslope_plus_1" dimension
19    INTEGER,SAVE :: timeindex ! current time index (for time-dependent fields)
20    INTEGER,PARAMETER :: length=100 ! size of tab_cntrl array
21   
22    INTERFACE get_field
23      MODULE PROCEDURE Get_field_r1,Get_field_r2,Get_field_r3
24    END INTERFACE get_field
25   
26    INTERFACE get_var
27      MODULE PROCEDURE get_var_r0,Get_var_r1,Get_var_r2,Get_var_r3
28    END INTERFACE get_var
29
30    INTERFACE put_field
31      MODULE PROCEDURE put_field_r1,put_field_r2,put_field_r3
32    END INTERFACE put_field
33
34    INTERFACE put_var
35      MODULE PROCEDURE put_var_r0,put_var_r1,put_var_r2,put_var_r3
36    END INTERFACE put_var
37
38    PUBLIC nid_start, length
39    PUBLIC get_field,get_var,put_field,put_var
40    PUBLIC inquire_dimension, inquire_dimension_length
41    PUBLIC inquire_field, inquire_field_ndims
42    PUBLIC open_startphy,close_startphy,open_restartphy,close_restartphy
43   
44CONTAINS
45
46  SUBROUTINE open_startphy(filename)
47  USE netcdf, only: NF90_OPEN, NF90_NOERR, NF90_NOWRITE, nf90_strerror
48  USE mod_phys_lmdz_para, only: is_master, bcast
49  IMPLICIT NONE
50    CHARACTER(LEN=*) :: filename
51    INTEGER          :: ierr
52
53    IF (is_master) THEN
54      ierr = NF90_OPEN (filename, NF90_NOWRITE, nid_start)
55      IF (ierr.NE.NF90_NOERR) THEN
56        write(*,*)'open_startphy: problem opening file '//trim(filename)
57        write(*,*)trim(nf90_strerror(ierr))
58        CALL abort_physic("open_startphy","Cannot open file",1)
59      ENDIF
60    ENDIF
61   
62    CALL bcast(nid_start) ! tell all procs about nid_start
63 
64  END SUBROUTINE open_startphy
65
66  SUBROUTINE close_startphy
67  USE netcdf, only: NF90_CLOSE
68  USE mod_phys_lmdz_para, only: is_master
69  IMPLICIT NONE
70    INTEGER          :: ierr
71
72    IF (is_master) THEN
73        ierr = NF90_CLOSE (nid_start)
74    ENDIF
75
76  END SUBROUTINE close_startphy
77
78
79  FUNCTION inquire_field(Field_name)
80  ! check if a given field is present in the input file
81  USE netcdf, only: NF90_INQ_VARID, NF90_NOERR
82  USE mod_phys_lmdz_para, only: is_master, bcast
83  IMPLICIT NONE
84    CHARACTER(LEN=*),INTENT(IN) :: Field_name
85    LOGICAL :: inquire_field
86    INTEGER :: varid
87    INTEGER :: ierr
88   
89    IF (is_master) THEN
90      ierr=NF90_INQ_VARID(nid_start,Field_name,varid)
91      IF (ierr==NF90_NOERR) THEN
92        Inquire_field=.TRUE.
93      ELSE
94        Inquire_field=.FALSE.
95      ENDIF
96    ENDIF
97
98    CALL bcast(inquire_field)
99
100  END FUNCTION inquire_field
101
102
103  FUNCTION inquire_field_ndims(Field_name)
104  ! give the number of dimensions of "Field_name" stored in the input file
105  USE netcdf, only: nf90_inq_varid, nf90_inquire_variable, &
106                    NF90_NOERR, nf90_strerror
107  USE mod_phys_lmdz_para, only: is_master, bcast
108  IMPLICIT NONE
109    CHARACTER(LEN=*),INTENT(IN) :: Field_name
110    INTEGER :: inquire_field_ndims
111    INTEGER :: varid
112    INTEGER :: ierr
113   
114    IF (is_master) THEN
115      ierr=nf90_inq_varid(nid_start,Field_name,varid)
116      ierr=nf90_inquire_variable(nid_start,varid,&
117                                  ndims=inquire_field_ndims)
118      IF (ierr.NE.NF90_NOERR) THEN
119        write(*,*)'inquire_field_ndims: problem obtaining ndims of '&
120                  //trim(field_name)
121        write(*,*)trim(nf90_strerror(ierr))
122        CALL abort_physic("inquire_field_ndims","Failed to get ndims",1)
123      ENDIF
124    ENDIF
125
126    CALL bcast(inquire_field_ndims)
127
128  END FUNCTION inquire_field_ndims
129
130
131  FUNCTION inquire_dimension(Field_name)
132  ! check if a given dimension is present in the input file
133  USE netcdf, only: nf90_inq_dimid, NF90_NOERR
134  USE mod_phys_lmdz_para, only: is_master, bcast
135  IMPLICIT NONE
136    CHARACTER(LEN=*),INTENT(IN) :: Field_name
137    LOGICAL :: inquire_dimension
138    INTEGER :: varid
139    INTEGER :: ierr
140   
141    IF (is_master) THEN
142      ierr=NF90_INQ_DIMID(nid_start,Field_name,varid)
143      IF (ierr==NF90_NOERR) THEN
144        Inquire_dimension=.TRUE.
145      ELSE
146        Inquire_dimension=.FALSE.
147      ENDIF
148    ENDIF
149
150    CALL bcast(inquire_dimension)
151
152  END FUNCTION inquire_dimension
153
154  FUNCTION inquire_dimension_length(Field_name)
155  ! give the length of the "Field_name" dimension stored in the input file
156  USE netcdf, only: nf90_inquire_dimension, nf90_inq_dimid, &
157                    NF90_NOERR, nf90_strerror
158  USE mod_phys_lmdz_para, only: is_master, bcast
159  IMPLICIT NONE
160    CHARACTER(LEN=*),INTENT(IN) :: Field_name
161    INTEGER :: inquire_dimension_length
162    INTEGER :: varid
163    INTEGER :: ierr
164   
165    IF (is_master) THEN
166      ierr=nf90_inq_dimid(nid_start,Field_name,varid)
167      ierr=nf90_inquire_dimension(nid_start,varid,&
168                                  len=inquire_dimension_length)
169      IF (ierr.NE.NF90_NOERR) THEN
170        write(*,*)'inquire_field_length: problem obtaining length of '&
171                  //trim(field_name)
172        write(*,*)trim(nf90_strerror(ierr))
173        CALL abort_physic("inquire_field_ndims","Failed to get length",1)
174      ENDIF
175    ENDIF
176
177    CALL bcast(inquire_dimension_length)
178
179  END FUNCTION inquire_dimension_length
180
181
182
183  SUBROUTINE Get_Field_r1(field_name,field,found,timeindex)
184  ! For a surface field
185  use mod_grid_phy_lmdz, only: klon_glo ! number of atmospheric columns (full grid)
186  IMPLICIT NONE
187    CHARACTER(LEN=*),INTENT(IN)    :: Field_name
188    REAL,INTENT(INOUT)               :: Field(:)
189    LOGICAL,INTENT(OUT),OPTIONAL   :: found
190    INTEGER,INTENT(IN),OPTIONAL    :: timeindex ! time index of sought data
191
192    integer :: edges(4), corners(4)
193
194    edges(1)=klon_glo
195    edges(2:4)=1
196    corners(1:4)=1
197    if (PRESENT(timeindex)) then
198      corners(2)=timeindex
199    endif
200
201    IF (PRESENT(found)) THEN
202      CALL Get_field_rgen(field_name,field,1,corners,edges,found)
203    ELSE
204      CALL Get_field_rgen(field_name,field,1,corners,edges)
205    ENDIF
206     
207  END SUBROUTINE Get_Field_r1
208 
209  SUBROUTINE Get_Field_r2(field_name,field,found,timeindex)
210  ! For a "3D" horizontal-vertical field
211  use mod_grid_phy_lmdz, only: klon_glo ! number of atmospheric columns (full grid)
212  IMPLICIT NONE
213    CHARACTER(LEN=*),INTENT(IN)    :: Field_name
214    REAL,INTENT(INOUT)               :: Field(:,:)
215    LOGICAL,INTENT(OUT),OPTIONAL   :: found
216    INTEGER,INTENT(IN),OPTIONAL    :: timeindex ! time index of sought data
217
218    integer :: edges(4), corners(4)
219
220    edges(1)=klon_glo
221    edges(2)=size(field,2)
222    edges(3:4)=1
223    corners(1:4)=1
224    if (PRESENT(timeindex)) then
225      corners(3)=timeindex
226    endif
227   
228    IF (PRESENT(found)) THEN
229      CALL Get_field_rgen(field_name,field,size(field,2),&
230                          corners,edges,found)
231    ELSE
232      CALL Get_field_rgen(field_name,field,size(field,2),&
233                          corners,edges)
234    ENDIF
235
236     
237  END SUBROUTINE Get_Field_r2
238 
239  SUBROUTINE Get_Field_r3(field_name,field,found,timeindex)
240  ! for a "4D" field surf/alt/??
241  use mod_grid_phy_lmdz, only: klon_glo ! number of atmospheric columns (full grid)
242  IMPLICIT NONE
243    CHARACTER(LEN=*),INTENT(IN)    :: Field_name
244    REAL,INTENT(INOUT)               :: Field(:,:,:)
245    LOGICAL,INTENT(OUT),OPTIONAL   :: found
246    INTEGER,INTENT(IN),OPTIONAL    :: timeindex ! time index of sought data
247
248    integer :: edges(4), corners(4)
249
250    edges(1)=klon_glo
251    edges(2)=size(field,2)
252    edges(3)=size(field,3)
253    edges(4)=1
254    corners(1:4)=1
255    if (PRESENT(timeindex)) then
256      corners(4)=timeindex
257    endif
258   
259    IF (PRESENT(found)) THEN
260      CALL Get_field_rgen(field_name,field,size(field,2)*size(field,3),&
261                          corners,edges,found)
262    ELSE
263      CALL Get_field_rgen(field_name,field,size(field,2)*size(field,3),&
264                          corners,edges)
265    ENDIF
266     
267  END SUBROUTINE Get_Field_r3
268 
269  SUBROUTINE Get_field_rgen(field_name,field,field_size, &
270                            corners,edges,found)
271  USE netcdf
272  USE dimphy
273  USE mod_grid_phy_lmdz
274  USE mod_phys_lmdz_para
275  IMPLICIT NONE
276    CHARACTER(LEN=*) :: Field_name
277    INTEGER          :: field_size
278    REAL             :: field(klon,field_size)
279    INTEGER,INTENT(IN) :: corners(4)
280    INTEGER,INTENT(IN) :: edges(4)
281    LOGICAL,OPTIONAL :: found
282   
283    REAL    :: field_glo(klon_glo,field_size)
284    LOGICAL :: tmp_found
285    INTEGER :: varid
286    INTEGER :: ierr
287   
288    IF (is_master) THEN
289 
290      ierr=NF90_INQ_VARID(nid_start,Field_name,varid)
291     
292      IF (ierr==NF90_NOERR) THEN
293        CALL body(field_glo)
294        tmp_found=.TRUE.
295      ELSE
296        tmp_found=.FALSE.
297      ENDIF
298   
299    ENDIF
300   
301    CALL bcast(tmp_found)
302
303    IF (tmp_found) THEN
304      CALL scatter(field_glo,field)
305    ENDIF
306   
307    IF (PRESENT(found)) THEN
308      found=tmp_found
309    ELSE
310      IF (.NOT. tmp_found) THEN
311        PRINT*, 'get_field_rgen: Field <'//field_name//'> not found'
312        CALL abort_physic("get_field_rgen","Failed to get field",1)
313      ENDIF
314    ENDIF
315 
316   
317    CONTAINS
318     
319     SUBROUTINE body(field_glo)
320       REAL :: field_glo(klon_glo*field_size)
321         ierr=NF90_GET_VAR(nid_start,varid,field_glo,corners,edges)
322         IF (ierr/=NF90_NOERR) THEN
323           ! La variable exist dans le fichier mais la lecture a echouee.
324           PRINT*, 'get_field_rgen: Failed reading <'//field_name//'>'
325
326!           IF (field_name=='CLWCON' .OR. field_name=='RNEBCON' .OR. field_name=='RATQS') THEN
327!              ! Essaye de lire le variable sur surface uniqument, comme fait avant
328!              field_glo(:)=0.
329!              ierr=NF90_GET_VAR(nid_start,varid,field_glo(1:klon_glo))
330!              IF (ierr/=NF90_NOERR) THEN
331!                 PRINT*, 'phyetat0: Lecture echouee aussi en 2D pour <'//field_name//'>'
332!                 CALL abort
333!              ELSE
334!                 PRINT*, 'phyetat0: La variable <'//field_name//'> lu sur surface seulement'!, selon ancien format, le reste mis a zero'
335!              END IF
336!           ELSE
337              CALL abort_physic("get_field_rgen","Failed to read field",1)
338!           ENDIF
339         ENDIF
340
341     END SUBROUTINE body
342
343  END SUBROUTINE Get_field_rgen
344
345
346  SUBROUTINE get_var_r0(var_name,var,found)
347  ! Get a scalar from input file
348  IMPLICIT NONE 
349    CHARACTER(LEN=*),INTENT(IN)  :: var_name
350    REAL,INTENT(INOUT)             :: var
351    LOGICAL,OPTIONAL,INTENT(OUT) :: found
352
353    REAL                         :: varout(1)
354   
355    IF (PRESENT(found)) THEN
356      CALL Get_var_rgen(var_name,varout,size(varout),found)
357    ELSE
358      CALL Get_var_rgen(var_name,varout,size(varout))
359    ENDIF
360    var=varout(1)
361 
362  END SUBROUTINE get_var_r0
363
364  SUBROUTINE get_var_r1(var_name,var,found)
365  ! Get a vector from input file
366  IMPLICIT NONE 
367    CHARACTER(LEN=*),INTENT(IN)  :: var_name
368    REAL,INTENT(INOUT)             :: var(:)
369    LOGICAL,OPTIONAL,INTENT(OUT) :: found
370   
371    IF (PRESENT(found)) THEN
372      CALL Get_var_rgen(var_name,var,size(var),found)
373    ELSE
374      CALL Get_var_rgen(var_name,var,size(var))
375    ENDIF
376 
377  END SUBROUTINE get_var_r1
378
379  SUBROUTINE get_var_r2(var_name,var,found)
380  ! Get a 2D field from input file
381  IMPLICIT NONE 
382    CHARACTER(LEN=*),INTENT(IN)  :: var_name
383    REAL,INTENT(OUT)             :: var(:,:)
384    LOGICAL,OPTIONAL,INTENT(OUT) :: found
385   
386    IF (PRESENT(found)) THEN
387      CALL Get_var_rgen(var_name,var,size(var),found)
388    ELSE
389      CALL Get_var_rgen(var_name,var,size(var))
390    ENDIF
391 
392  END SUBROUTINE get_var_r2
393
394  SUBROUTINE get_var_r3(var_name,var,found)
395  ! Get a 3D field frominput file
396  IMPLICIT NONE 
397    CHARACTER(LEN=*),INTENT(IN)  :: var_name
398    REAL,INTENT(INOUT)             :: var(:,:,:)
399    LOGICAL,OPTIONAL,INTENT(OUT) :: found
400   
401    IF (PRESENT(found)) THEN
402      CALL Get_var_rgen(var_name,var,size(var),found)
403    ELSE
404      CALL Get_var_rgen(var_name,var,size(var))
405    ENDIF
406 
407  END SUBROUTINE get_var_r3
408
409  SUBROUTINE Get_var_rgen(var_name,var,var_size,found)
410  USE netcdf
411  USE dimphy
412  USE mod_grid_phy_lmdz
413  USE mod_phys_lmdz_para
414  IMPLICIT NONE
415    CHARACTER(LEN=*) :: var_name
416    INTEGER          :: var_size
417    REAL             :: var(var_size)
418    LOGICAL,OPTIONAL :: found
419   
420    LOGICAL :: tmp_found
421    INTEGER :: varid
422    INTEGER :: ierr
423   
424    IF (is_mpi_root .AND. is_omp_root) THEN
425 
426      ierr=NF90_INQ_VARID(nid_start,var_name,varid)
427      IF (ierr==NF90_NOERR) THEN
428        ierr=NF90_GET_VAR(nid_start,varid,var)
429        IF (ierr/=NF90_NOERR) THEN
430          PRINT*, 'phyetat0: Failed loading <'//trim(var_name)//'>'
431          CALL abort_physic("get_var_rgen","Failed to read variable",1)
432        ENDIF
433        tmp_found=.TRUE.
434      ELSE
435        tmp_found=.FALSE.
436      ENDIF
437   
438    ENDIF
439   
440    CALL bcast(tmp_found)
441
442    IF (tmp_found) THEN
443      CALL bcast(var)
444    ENDIF
445   
446    IF (PRESENT(found)) THEN
447      found=tmp_found
448    ELSE
449      IF (.NOT. tmp_found) THEN
450        PRINT*, 'phyetat0: Variable <'//trim(var_name)//'> not found'
451        CALL abort_physic("get_var_rgen","Failed to read variable",1)
452      ENDIF
453    ENDIF
454
455  END SUBROUTINE Get_var_rgen
456
457!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
458
459  SUBROUTINE open_restartphy(filename)
460  USE netcdf, only: NF90_CREATE, NF90_CLOBBER, NF90_64BIT_OFFSET, &
461                    NF90_NOERR, nf90_strerror, &
462                    NF90_PUT_ATT, NF90_GLOBAL, NF90_DEF_DIM, &
463                    NF90_UNLIMITED, NF90_ENDDEF, &
464                    NF90_WRITE, NF90_OPEN
465  USE mod_phys_lmdz_para, only: is_master
466  USE mod_grid_phy_lmdz, only: klon_glo
467  USE dimphy, only: klev, klevp1
468  USE tracer_mod, only: nqmx
469  USE comsoil_h, only: nsoilmx
470  USE comslope_mod, only: nslope
471  IMPLICIT NONE
472    CHARACTER(LEN=*),INTENT(IN) :: filename
473    INTEGER                     :: ierr
474    LOGICAL,SAVE :: already_created=.false.
475    LOGICAL,SAVE :: already_created_restarfi_1d=.false.
476
477      if(klon_glo.eq.1 .and. filename.eq."restartfi.nc" .and. .not.already_created_restarfi_1d) then
478        already_created=.false.
479        already_created_restarfi_1d=.true.
480        first_time1d=.true.
481      endif
482   
483    IF (is_master) THEN
484      if (.not.already_created) then
485        ! At the very first call, create the file
486        ierr=NF90_CREATE(filename,IOR(NF90_CLOBBER,NF90_64BIT_OFFSET), &
487                          nid_restart)
488        IF (ierr/=NF90_NOERR) THEN
489          write(*,*)'open_restartphy: problem creating file '//trim(filename)
490          write(*,*)trim(nf90_strerror(ierr))
491          CALL abort_physic("open_restartphy","Failed creating file",1)
492        ENDIF
493        already_created=.true.
494      else
495        ! Just open the file
496        ierr=NF90_OPEN(filename,NF90_WRITE,nid_restart)
497        IF (ierr/=NF90_NOERR) THEN
498          write(*,*)'open_restartphy: problem opening file '//trim(filename)
499          write(*,*)trim(nf90_strerror(ierr))
500          CALL abort_physic("open_restartphy","Failed opening file",1)
501        ENDIF
502        return
503      endif ! of if (.not.already_created)
504
505      ierr=NF90_PUT_ATT(nid_restart,NF90_GLOBAL,"title",&
506                        "Physics start file")
507      IF (ierr/=NF90_NOERR) THEN
508        write(*,*)'open_restartphy: problem writing title '
509        write(*,*)trim(nf90_strerror(ierr))
510      ENDIF
511
512      ierr=NF90_DEF_DIM(nid_restart,"index",length,idim1)
513      IF (ierr/=NF90_NOERR) THEN
514        write(*,*)'open_restartphy: problem defining index dimension '
515        write(*,*)trim(nf90_strerror(ierr))
516        CALL abort_physic("open_restartphy","Failed defining index",1)
517      ENDIF
518     
519      ierr=NF90_DEF_DIM(nid_restart,"physical_points",klon_glo,idim2)
520      IF (ierr/=NF90_NOERR) THEN
521        write(*,*)'open_restartphy: problem defining physical_points dimension '
522        write(*,*)trim(nf90_strerror(ierr))
523        CALL abort_physic("open_restartphy","Failed defining physical_points",1)
524      ENDIF
525     
526      ierr=NF90_DEF_DIM(nid_restart,"subsurface_layers",nsoilmx,idim3)
527      IF (ierr/=NF90_NOERR) THEN
528        write(*,*)'open_restartphy: problem defining subsurface_layers dimension '
529        write(*,*)trim(nf90_strerror(ierr))
530        CALL abort_physic("open_restartphy","Failed defining subsurface_layers",1)
531      ENDIF
532     
533      ierr=NF90_DEF_DIM(nid_restart,"nlayer_plus_1",klevp1,idim4)
534      IF (ierr/=NF90_NOERR) THEN
535        write(*,*)'open_restartphy: problem defining nlayer_plus_1 dimension '
536        write(*,*)trim(nf90_strerror(ierr))
537        CALL abort_physic("open_restartphy","Failed defining nlayer_plus_1",1)
538      ENDIF
539     
540      if (nqmx.ne.0) then
541        ierr=NF90_DEF_DIM(nid_restart,"number_of_advected_fields",nqmx,idim5)
542      else
543        ! pretend nqmx=1 because 0 size implies unlimited dimension for NetCDF
544        ierr=NF90_DEF_DIM(nid_restart,"number_of_advected_fields",1,idim5)
545      endif
546      IF (ierr/=NF90_NOERR) THEN
547        write(*,*)'open_restartphy: problem defining number_of_advected_fields dimension '
548        write(*,*)trim(nf90_strerror(ierr))
549        CALL abort_physic("open_restartphy","Failed defining number_of_advected_fields",1)
550      ENDIF
551
552      ierr=NF90_DEF_DIM(nid_restart,"nlayer",klev,idim6)
553      IF (ierr/=NF90_NOERR) THEN
554        write(*,*)'open_restartphy: problem defining nlayer dimension '
555        write(*,*)trim(nf90_strerror(ierr))
556        CALL abort_physic("open_restartphy","Failed defining nlayer",1)
557      ENDIF
558     
559      ierr=NF90_DEF_DIM(nid_restart,"Time",NF90_UNLIMITED,idim7)
560      IF (ierr/=NF90_NOERR) THEN
561        write(*,*)'open_restartphy: problem defining Time dimension '
562        write(*,*)trim(nf90_strerror(ierr))
563        CALL abort_physic("open_restartphy","Failed defining Time",1)
564      ENDIF
565
566      ierr=NF90_DEF_DIM(nid_restart,"nslope",nslope,idim8)
567      IF (ierr/=NF90_NOERR) THEN
568        write(*,*)'phyredem: problem defining nslope dimension'
569        write(*,*)trim(nf90_strerror(ierr))
570        CALL ABORT
571      ENDIF
572
573      ierr=NF90_DEF_DIM(nid_restart,"inter_slope",nslope+1,idim9)
574      IF (ierr/=NF90_NOERR) THEN
575        write(*,*)'phyredem: problem defining inter slope dimension'
576        write(*,*)trim(nf90_strerror(ierr))
577        CALL ABORT
578      ENDIF
579
580      ierr=NF90_ENDDEF(nid_restart)
581      IF (ierr/=NF90_NOERR) THEN
582        write(*,*)'open_restartphy: problem ending definition mode '
583        write(*,*)trim(nf90_strerror(ierr))
584        CALL abort_physic("open_restartphy","Failed ending definition mode",1)
585      ENDIF
586    ENDIF
587
588  END SUBROUTINE open_restartphy
589
590  SUBROUTINE close_restartphy
591  USE netcdf, only: NF90_CLOSE
592  USE mod_phys_lmdz_para, only: is_master
593  IMPLICIT NONE
594    INTEGER          :: ierr
595
596    IF (is_master) THEN
597      ierr = NF90_CLOSE (nid_restart)
598    ENDIF
599 
600  END SUBROUTINE close_restartphy
601
602  SUBROUTINE put_field_r1(field_name,title,field,time)
603  ! For a surface field
604  IMPLICIT NONE
605  CHARACTER(LEN=*),INTENT(IN)    :: field_name
606  CHARACTER(LEN=*),INTENT(IN)    :: title
607  REAL,INTENT(IN)                :: field(:)
608  REAL,OPTIONAL,INTENT(IN)       :: time
609 
610  IF (present(time)) THEN
611    ! if timeindex is present, it is a time-dependent variable
612    CALL put_field_rgen(field_name,title,field,1,time)
613  ELSE
614    CALL put_field_rgen(field_name,title,field,1)
615  ENDIF
616 
617  END SUBROUTINE put_field_r1
618
619  SUBROUTINE put_field_r2(field_name,title,field,time)
620  ! For a "3D" horizontal-vertical field
621  IMPLICIT NONE
622  CHARACTER(LEN=*),INTENT(IN)    :: field_name
623  CHARACTER(LEN=*),INTENT(IN)    :: title
624  REAL,INTENT(IN)                :: field(:,:)
625  REAL,OPTIONAL,INTENT(IN)       :: time
626 
627  IF (present(time)) THEN
628    ! if timeindex is present, it is a time-dependent variable
629    CALL put_field_rgen(field_name,title,field,size(field,2),time)
630  ELSE
631    CALL put_field_rgen(field_name,title,field,size(field,2))
632  ENDIF
633 
634  END SUBROUTINE put_field_r2
635
636  SUBROUTINE put_field_r3(field_name,title,field,time)
637  ! For a "4D" field surf/alt/??
638  IMPLICIT NONE
639  CHARACTER(LEN=*),INTENT(IN)    :: field_name
640  CHARACTER(LEN=*),INTENT(IN)    :: title
641  REAL,INTENT(IN)                :: field(:,:,:)
642  REAL,OPTIONAL,INTENT(IN)       :: time
643 
644  IF (present(time)) THEN
645    ! if timeindex is present, it is a time-dependent variable
646    CALL put_field_rgen(field_name,title,field,size(field,2)*size(field,3),&
647                        time)
648  ELSE 
649    CALL put_field_rgen(field_name,title,field,size(field,2)*size(field,3))
650  ENDIF
651 
652  END SUBROUTINE put_field_r3
653 
654  SUBROUTINE put_field_rgen(field_name,title,field,field_size,time)
655  USE netcdf
656  USE dimphy
657  USE comsoil_h, only: nsoilmx
658  USE comslope_mod, ONLY: nslope
659  USE mod_grid_phy_lmdz
660  USE mod_phys_lmdz_para
661  IMPLICIT NONE
662  CHARACTER(LEN=*),INTENT(IN)    :: field_name
663  CHARACTER(LEN=*),INTENT(IN)    :: title
664  INTEGER,INTENT(IN)             :: field_size
665  REAL,INTENT(IN)                :: field(klon,field_size)
666  REAL,OPTIONAL,INTENT(IN)       :: time
667 
668  REAL                           :: field_glo(klon_glo,field_size)
669  REAL                           :: field_glo_reshape(klon_glo,nsoilmx,nslope,timeindex)
670  INTEGER                        :: ierr
671  INTEGER                        :: nvarid
672  INTEGER                        :: idim
673   
674    CALL gather(field,field_glo)
675   
676    IF (is_master) THEN
677
678      IF (field_size==1) THEN
679        ! input is a 1D "surface field" array
680        if (.not.present(time)) then ! for a time-independent field
681          ierr=NF90_REDEF(nid_restart)
682#ifdef NC_DOUBLE
683          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
684                            (/idim2/),nvarid)
685#else
686          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
687                            (/idim2/),nvarid)
688#endif
689          if (ierr.ne.NF90_NOERR) then
690            write(*,*)"put_field_rgen error: failed to define "//trim(field_name)
691            write(*,*)trim(nf90_strerror(ierr))
692          endif
693          IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
694          ierr=NF90_ENDDEF(nid_restart)
695          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo)
696        else
697          ! check if the variable has already been defined:
698          ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid)
699          if (ierr/=NF90_NOERR) then ! variable not found, define it
700            ierr=NF90_REDEF(nid_restart)
701#ifdef NC_DOUBLE
702            if(field_name.eq. "tauscaling" .or. field_name.eq. "totcloudfrac" .or. field_name.eq. "wstar" ) then
703            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
704                              (/idim2,idim7/),nvarid)
705            else
706            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
707                              (/idim2,idim8,idim7/),nvarid)
708            endif
709#else
710            if(field_name.eq. "tauscaling" .or. field_name.eq. "totcloudfrac" .or. field_name.eq. "wstar" ) then
711            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
712                              (/idim2,idim7/),nvarid)
713            else
714            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
715                              (/idim2,idim8,idim7/),nvarid)
716            endif
717#endif
718            if (ierr.ne.NF90_NOERR) then
719              write(*,*)"put_field_rgen error: failed to define "//trim(field_name)
720              write(*,*)trim(nf90_strerror(ierr))
721            endif
722            IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
723            ierr=NF90_ENDDEF(nid_restart)
724          endif
725          ! Write the variable
726          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,&
727                            start=(/1,timeindex/))
728        endif ! of if (.not.present(timeindex))
729
730      ELSE IF (field_size==klev) THEN
731        ! input is a 2D "atmospheric field" array
732        if (.not.present(time)) then ! for a time-independent field
733          ierr=NF90_REDEF(nid_restart)
734#ifdef NC_DOUBLE
735          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
736                            (/idim2,idim6/),nvarid)
737#else
738          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
739                            (/idim2,idim6/),nvarid)
740#endif
741          if (ierr.ne.NF90_NOERR) then
742            write(*,*)"put_field_rgen error: failed to define "//trim(field_name)
743            write(*,*)trim(nf90_strerror(ierr))
744          endif
745          IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
746          ierr=NF90_ENDDEF(nid_restart)
747          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo)
748        else
749          ! check if the variable has already been defined:
750          ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid)
751          if (ierr/=NF90_NOERR) then ! variable not found, define it
752            ierr=NF90_REDEF(nid_restart)
753#ifdef NC_DOUBLE
754            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
755                              (/idim2,idim6,idim7/),nvarid)
756#else
757            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
758                              (/idim2,idim6,idim7/),nvarid)
759#endif
760            if (ierr.ne.NF90_NOERR) then
761              write(*,*)"put_field_rgen error: failed to define "//trim(field_name)
762              write(*,*)trim(nf90_strerror(ierr))
763            endif
764            IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
765            ierr=NF90_ENDDEF(nid_restart)
766          endif
767          ! Write the variable
768          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,&
769                            start=(/1,1,timeindex/))
770        endif ! of if (.not.present(time))
771
772      ELSE IF (field_size==klevp1) THEN
773        ! input is a 2D "interlayer atmospheric field" array
774        if (.not.present(time)) then ! for a time-independent field
775          ierr=NF90_REDEF(nid_restart)
776#ifdef NC_DOUBLE
777          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
778                            (/idim2,idim4/),nvarid)
779#else
780          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
781                            (/idim2,idim4/),nvarid)
782#endif
783          if (ierr.ne.NF90_NOERR) then
784            write(*,*)"put_field_rgen error: failed to define "//trim(field_name)
785            write(*,*)trim(nf90_strerror(ierr))
786          endif
787          IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
788          ierr = NF90_ENDDEF(nid_restart)
789          ierr = NF90_PUT_VAR(nid_restart,nvarid,field_glo)
790        else
791          ! check if the variable has already been defined:
792          ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid)
793          if (ierr/=NF90_NOERR) then ! variable not found, define it
794            ierr=NF90_REDEF(nid_restart)
795#ifdef NC_DOUBLE
796            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
797                              (/idim2,idim4,idim7/),nvarid)
798#else
799            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
800                              (/idim2,idim4,idim7/),nvarid)
801#endif
802            if (ierr.ne.NF90_NOERR) then
803              write(*,*)"put_field_rgen error: failed to define "//trim(field_name)
804              write(*,*)trim(nf90_strerror(ierr))
805            endif
806            IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
807            ierr=NF90_ENDDEF(nid_restart)
808          endif
809          ! Write the variable
810          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,&
811                            start=(/1,1,timeindex/))
812        endif ! of if (.not.present(timeindex))
813
814      ELSE IF (field_size==nsoilmx .or. field_size==nsoilmx*nslope) THEN
815        ! input is a 2D "subsurface field" array
816        if (.not.present(time)) then ! for a time-independent field
817          ierr = NF90_REDEF(nid_restart)
818#ifdef NC_DOUBLE
819          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
820                            (/idim2,idim3/),nvarid)
821#else
822          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
823                            (/idim2,idim3/),nvarid)
824#endif
825          if (ierr.ne.NF90_NOERR) then
826            write(*,*)"put_field_rgen error: failed to define "//trim(field_name)
827            write(*,*)trim(nf90_strerror(ierr))
828          endif
829          IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
830          ierr = NF90_ENDDEF(nid_restart)
831          ierr = NF90_PUT_VAR(nid_restart,nvarid,field_glo)
832        else
833          ! check if the variable has already been defined:
834          ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid)
835          if (ierr/=NF90_NOERR) then ! variable not found, define it
836            ierr=NF90_REDEF(nid_restart)
837#ifdef NC_DOUBLE
838            if(field_name.eq. "tsoil".or. field_name.eq. "inertiesoil") then
839            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
840                              (/idim2,idim3,idim8,idim7/),nvarid)
841            else
842            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
843                              (/idim2,idim3,idim7/),nvarid)
844            endif
845#else
846            if(field_name.eq. "tsoil".or. field_name.eq. "inertiesoil") then
847            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
848                              (/idim2,idim3,idim8,idim7/),nvarid)
849            else
850            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
851                              (/idim2,idim3,idim7/),nvarid)
852            endif
853#endif
854           if (ierr.ne.NF90_NOERR) then
855              write(*,*)"put_field_rgen error: failed to define "//trim(field_name)
856              write(*,*)trim(nf90_strerror(ierr))
857            endif
858            IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
859            ierr=NF90_ENDDEF(nid_restart)
860          endif
861          ! Write the variable
862          if(field_name.eq. "tsoil" .or. field_name.eq. "inertiesoil") then
863
864          field_glo_reshape=RESHAPE(field_glo, (/klon_glo, nsoilmx, nslope, timeindex/))
865
866          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo_reshape,&
867                            start=(/1,1,1,timeindex/))
868          else
869          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,&
870                            start=(/1,1,timeindex/))
871          endif
872        endif ! of if (.not.present(time))
873
874      ELSE IF (field_size==nslope) THEN
875        ! input is a 2D "subsurface field" array
876        if (.not.present(time)) then ! for a time-independent field
877          ierr = NF90_REDEF(nid_restart)
878#ifdef NC_DOUBLE
879          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
880                            (/idim2,idim8/),nvarid)
881#else
882          ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
883                            (/idim2,idim8/),nvarid)
884#endif
885          if (ierr.ne.NF90_NOERR) then
886            write(*,*)"put_field_rgen error: failed to define"//trim(field_name)
887            write(*,*)trim(nf90_strerror(ierr))
888          endif
889          IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
890          ierr = NF90_ENDDEF(nid_restart)
891          ierr = NF90_PUT_VAR(nid_restart,nvarid,field_glo)
892        else
893          ! check if the variable has already been defined:
894          ierr=NF90_INQ_VARID(nid_restart,field_name,nvarid)
895          if (ierr/=NF90_NOERR) then ! variable not found, define it
896            ierr=NF90_REDEF(nid_restart)
897#ifdef NC_DOUBLE
898            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
899                              (/idim2,idim8,idim7/),nvarid)
900#else
901            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
902                              (/idim2,idim8,idim7/),nvarid)
903#endif
904           if (ierr.ne.NF90_NOERR) then
905              write(*,*)"put_field_rgen error: failed to define"//trim(field_name)
906              write(*,*)trim(nf90_strerror(ierr))
907            endif
908            IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
909            ierr=NF90_ENDDEF(nid_restart)
910          endif
911          ! Write the variable
912          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,&
913                            start=(/1,1,timeindex/))
914
915        endif ! of if (.not.present(time))
916
917      ELSE
918        PRINT *, "Error phyredem(put_field_rgen) : wrong dimension for ",trim(field_name)
919        write(*,*) "  field_size =",field_size
920        CALL abort_physic("put_field_rgen","wrong field dimension",1)
921      ENDIF
922
923      ! Check the writting of field to file went OK
924      if (ierr.ne.NF90_NOERR) then
925        write(*,*) " Error phyredem(put_field_rgen) : failed writing ",trim(field_name)
926        write(*,*)trim(nf90_strerror(ierr))
927        CALL abort_physic("put_field_rgen","Failed writing field",1)
928      endif
929
930    ENDIF ! of IF (is_master)
931   
932  END SUBROUTINE put_field_rgen 
933 
934  SUBROUTINE put_var_r0(var_name,title,var)
935  ! Put a scalar in file
936   IMPLICIT NONE
937     CHARACTER(LEN=*),INTENT(IN) :: var_name
938     CHARACTER(LEN=*),INTENT(IN) :: title
939     REAL,INTENT(IN)             :: var
940     REAL                        :: varin(1)
941     
942     varin(1)=var
943     
944     CALL put_var_rgen(var_name,title,varin,size(varin))
945
946  END SUBROUTINE put_var_r0
947
948
949  SUBROUTINE put_var_r1(var_name,title,var)
950  ! Put a vector in file
951   IMPLICIT NONE
952     CHARACTER(LEN=*),INTENT(IN) :: var_name
953     CHARACTER(LEN=*),INTENT(IN) :: title
954     REAL,INTENT(IN)             :: var(:)
955     
956     CALL put_var_rgen(var_name,title,var,size(var))
957
958  END SUBROUTINE put_var_r1
959 
960  SUBROUTINE put_var_r2(var_name,title,var)
961  ! Put a 2D field in file
962   IMPLICIT NONE
963     CHARACTER(LEN=*),INTENT(IN) :: var_name
964     CHARACTER(LEN=*),INTENT(IN) :: title
965     REAL,INTENT(IN)             :: var(:,:)
966     
967     CALL put_var_rgen(var_name,title,var,size(var))
968
969  END SUBROUTINE put_var_r2     
970 
971  SUBROUTINE put_var_r3(var_name,title,var)
972  ! Put a 3D field in file
973   IMPLICIT NONE
974     CHARACTER(LEN=*),INTENT(IN) :: var_name
975     CHARACTER(LEN=*),INTENT(IN) :: title
976     REAL,INTENT(IN)             :: var(:,:,:)
977     
978     CALL put_var_rgen(var_name,title,var,size(var))
979
980  END SUBROUTINE put_var_r3
981
982  SUBROUTINE put_var_rgen(var_name,title,var,var_size)
983  USE netcdf, only: NF90_REDEF, NF90_DEF_VAR, NF90_ENDDEF, NF90_PUT_VAR, &
984                    NF90_FLOAT, NF90_DOUBLE, &
985                    NF90_PUT_ATT, NF90_NOERR, nf90_strerror, &
986                    nf90_inq_dimid, nf90_inquire_dimension, NF90_INQ_VARID
987  USE comsoil_h, only: nsoilmx
988  USE comslope_mod, only: nslope
989  USE mod_phys_lmdz_para, only: is_master
990  IMPLICIT NONE
991     CHARACTER(LEN=*),INTENT(IN) :: var_name
992     CHARACTER(LEN=*),INTENT(IN) :: title
993     INTEGER,INTENT(IN)          :: var_size
994     REAL,INTENT(IN)             :: var(var_size)
995     
996     INTEGER :: ierr
997     INTEGER :: nvarid
998     INTEGER :: idim1d
999     logical,save :: firsttime=.true.
1000         
1001    IF (is_master) THEN
1002
1003      IF (var_name=="Time") THEN
1004        ! Very specific case of "Time" variable
1005        if (firsttime .or. first_time1d) then
1006          ! Create the "Time variable"
1007          ierr=NF90_REDEF(nid_restart)
1008#ifdef NC_DOUBLE
1009          ierr=NF90_DEF_VAR(nid_restart,var_name,NF90_DOUBLE,&
1010                            (/idim7/),nvarid)
1011#else
1012          ierr=NF90_DEF_VAR(nid_restart,var_name,NF90_FLOAT,&
1013                            (/idim7/),nvarid)
1014#endif
1015          IF (LEN_TRIM(title) > 0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
1016          ierr=NF90_ENDDEF(nid_restart)
1017          firsttime=.false.
1018          if(first_time1d) first_time1d=.false.
1019        endif
1020        ! Append "time" value
1021        ! get current length of "Time" dimension
1022        ierr=nf90_inq_dimid(nid_restart,var_name,idim1d)
1023        ierr=nf90_inquire_dimension(nid_restart,idim1d,len=timeindex)
1024        timeindex=timeindex+1
1025        ierr=NF90_INQ_VARID(nid_restart,var_name,nvarid)
1026        ierr=NF90_PUT_VAR(nid_restart,nvarid,var,&
1027                           start=(/timeindex/))
1028        IF (ierr/=NF90_NOERR) THEN
1029          write(*,*)'put_var_rgen: problem writing Time'
1030          write(*,*)trim(nf90_strerror(ierr))
1031          CALL abort_physic("get_var_rgen","Failed to write Time",1)
1032        ENDIF
1033        return ! nothing left to do
1034      ELSEIF (var_size==length) THEN
1035        ! We know it is a "controle" kind of 1D array
1036        idim1d=idim1
1037      ELSEIF (var_size==nsoilmx) THEN
1038        ! We know it is an  "mlayer" kind of 1D array
1039        idim1d=idim3
1040      ELSEIF (var_size==nslope+1) THEN
1041        ! We know it is an  "inter slope" kind of 1D array
1042        idim1d=idim9
1043      ELSE
1044        PRINT *, "put_var_rgen error : wrong dimension"
1045        write(*,*) "  var_size =",var_size
1046        CALL abort_physic("get_var_rgen","Wrong variable dimension",1)
1047
1048      ENDIF ! of IF (var_size==length) THEN
1049
1050      ! Swich to NetCDF define mode
1051      ierr=NF90_REDEF (nid_restart)
1052      ! Define the variable
1053#ifdef NC_DOUBLE
1054      ierr=NF90_DEF_VAR(nid_restart,var_name,NF90_DOUBLE,(/idim1d/),nvarid)
1055#else
1056      ierr=NF90_DEF_VAR(nid_restart,var_name,NF90_FLOAT,(/idim1d/),nvarid)
1057#endif
1058     
1059      ! Add a "title" attribute
1060      IF (LEN_TRIM(title)>0) ierr=NF90_PUT_ATT(nid_restart,nvarid,"title",title)
1061      ! Swich out of define mode
1062      ierr=NF90_ENDDEF(nid_restart)
1063      ! Write variable to file
1064      ierr=NF90_PUT_VAR(nid_restart,nvarid,var)
1065      IF (ierr/=NF90_NOERR) THEN
1066        write(*,*)'put_var_rgen: problem writing '//trim(var_name)
1067        write(*,*)trim(nf90_strerror(ierr))
1068        CALL abort_physic("put_var_rgen","Failed writing variable",1)
1069      ENDIF
1070    ENDIF ! of IF (is_master)
1071   
1072  END SUBROUTINE put_var_rgen     
1073
1074END MODULE iostart
Note: See TracBrowser for help on using the repository browser.