source: trunk/LMDZ.PLUTO/libf/phypluto/writediagspecVI.F @ 3955

Last change on this file since 3955 was 3749, checked in by afalco, 8 months ago

Pluto: ecritphy changed into diagfi_output_rate (followup of generic model), which defines the output rate of the diagfi in physical timesteps rather than dynamical ones.
AF

File size: 11.4 KB
RevLine 
[3184]1      subroutine writediagspecVI(ngrid,nom,titre,unite,dimpx,px)
2
[3749]3!  Ecriture de variables diagnostiques au choix dans la physique
[3184]4!  dans un fichier NetCDF nomme  'diagfi'. Ces variables peuvent etre
5!  3d (ex : temperature), 2d (ex : temperature de surface), ou
6!  0d (pour un scalaire qui ne depend que du temps : ex : la longitude
7!  solaire)
[3749]8!  Dans la version 2000, la periode d'ecriture est celle de
9!  "diagfi_output_rate " regle dans le fichier de controle de run :  run.def
[3184]10!
11!    writediagfi peut etre appele de n'importe quelle subroutine
12!    de la physique, plusieurs fois. L'initialisation et la creation du
13!    fichier se fait au tout premier appel.
14!
15! WARNING : les variables dynamique (u,v,t,q,ps)
16!  sauvees par writediagfi avec une
[3749]17! date donnee sont legerement differentes que dans le fichier histoire car
[3184]18! on ne leur a pas encore ajoute de la dissipation et de la physique !!!
19! IL est  RECOMMANDE d'ajouter les tendance physique a ces variables
20! avant l'ecriture dans diagfi (cf. physiq.F)
21!
[3749]22!
[3184]23!  parametres (input) :
24!  ----------
25!      ngrid : nombres de point ou est calcule la physique
26!                (ngrid = 2+(jjm-1)*iim - 1/jjm)
27!                 (= nlon ou klon dans la physique terrestre)
[3749]28!
[3184]29!      unit : unite logique du fichier de sortie (toujours la meme)
30!      nom  : nom de la variable a sortir (chaine de caracteres)
31!      titre: titre de la variable (chaine de caracteres)
32!      unite : unite de la variable (chaine de caracteres)
33!      px : variable a sortir (real 0, 2, ou 3d)
34!      dimpx : dimension de px : 0, 2, ou 3 dimensions
35!
36!=================================================================
37!
38!      This is a modified version that accepts spectrally varying input
39!      RW (2010)
40!
41!=================================================================
[3749]42
[3184]43! Addition by RW (2010) to allow OSR to be saved in .nc format
44      use radinc_h, only : L_NSPECTV
45      use geometry_mod, only: cell_area
46      use mod_phys_lmdz_para, only : is_mpi_root, is_master, gather
47      use mod_grid_phy_lmdz, only : klon_glo, Grid1Dto2D_glo,
48     &                              nbp_lon, nbp_lat, grid_type,
49     &                              unstructured
[3749]50      use time_phylmdz_mod, only: diagfi_output_rate,dtphys,daysec
51      use time_phylmdz_mod, only: day_ini
[3184]52      use callkeys_mod, only: iradia
53
54      implicit none
55
56      include "netcdf.inc"
57
58! Arguments on input:
59      integer ngrid
60      character (len=*) :: nom,titre,unite
61      integer dimpx
62      real px(ngrid,L_NSPECTV)
63
64! Local variables:
65
66!      real dx3(iip1,jjp1,llm) ! to store a 3D data set
67!      real dx2(iip1,jjp1)     ! to store a 2D (surface) data set
68!      real dx0
69
[3749]70      integer isample
[3184]71      integer ierr
72      integer iq
73      integer i,j,l,zmax , ig0
74
75      integer,save :: zitau=0
76      character(len=20),save :: firstnom='1234567890'
77      real,save :: date
78!$OMP THREADPRIVATE(firstnom,zitau,date)
79
80! Ajouts
81      integer, save :: ntime=0
82!$OMP THREADPRIVATE(ntime)
83      integer :: idim,varid
84      integer :: nid
85      character (len =50):: fichnom
86      integer, dimension(4) :: id
87      integer, dimension(4) :: edges,corner
88
89      real area((nbp_lon+1),nbp_lat)
90! added by RDW for OSR output
91      real dx3(nbp_lon+1,nbp_lat,L_NSPECTV) ! to store the data set
92      real dx3_1d(1,L_NSPECTV) ! to store the data with 1D model
93
94#ifdef CPP_PARA
95! Added to work in parallel mode
96      real dx3_glop(klon_glo,L_NSPECTV)
97      real dx3_glo(nbp_lon,nbp_lat,L_NSPECTV) ! to store a global 3D data set
98      real areafi_glo(klon_glo) ! mesh area on global physics grid
99#else
100      real areafi_glo(ngrid) ! mesh area on global physics grid
101#endif
[3749]102      if (grid_type == unstructured) then
[3184]103        return
104      endif
105
106!***************************************************************
107!Sortie des variables au rythme voulu
108
[3749]109      isample = diagfi_output_rate*iradia ! sortie au rythme de diagfi_output_rate
[3184]110!EM+JL if the spetra need to be output more frequently, need to define a ecritSpec...
[3749]111!     isample = diagfi_output_rate  ! sortie a tous les pas physique
[3184]112
113!***************************************************************
114
115! Initialisation of 'firstnom' and create/open the "diagfi.nc" NetCDF file
116! ------------------------------------------------------------------------
117! (Au tout premier appel de la subroutine durant le run.)
118
119      fichnom="diagspecVI.nc"
120
121      if (firstnom.eq.'1234567890') then ! .true. for the very first call
122      !  to this subroutine; now set 'firstnom'
123         firstnom = nom
124         ! just to be sure, check that firstnom is large enough to hold nom
125         if (len_trim(firstnom).lt.len_trim(nom)) then
126           write(*,*) "writediagfi: Error !!!"
127           write(*,*) "   firstnom string not long enough!!"
128           write(*,*) "   increase its size to at least ",len_trim(nom)
129           stop
130         endif
131
132#ifdef CPP_PARA
133          ! Gather cell_area() mesh area on physics grid
134          call Gather(cell_area,areafi_glo)
135#else
136         areafi_glo(:)=cell_area(:)
137#endif
138         ! Create the NetCDF file
139         if (is_master) then
140         ierr = NF_CREATE(fichnom, NF_CLOBBER, nid)
141         ! Define the 'Time' dimension
142         ierr = nf_def_dim(nid,"Time",NF_UNLIMITED,idim)
143         ! Define the 'Time' variable
144#ifdef NC_DOUBLE
145         ierr = NF_DEF_VAR (nid, "Time", NF_DOUBLE, 1, idim,varid)
146#else
147         ierr = NF_DEF_VAR (nid, "Time", NF_FLOAT, 1, idim,varid)
148#endif
149         ! Add a long_name attribute
150         ierr = NF_PUT_ATT_TEXT (nid, varid, "long_name",
151     .          4,"Time")
152         ! Add a units attribute
153         ierr = NF_PUT_ATT_TEXT(nid, varid,'units',29,
154     .          "days since 0000-00-0 00:00:00")
155         ! Switch out of NetCDF Define mode
156         ierr = NF_ENDDEF(nid)
157
158         ! Build area()
159         IF (klon_glo>1) THEN
160          do i=1,nbp_lon+1 ! poles
161           ! divide at the poles by nbp_lon
162           area(i,1)=areafi_glo(1)/nbp_lon
163           area(i,nbp_lat)=areafi_glo(klon_glo)/nbp_lon
164          enddo
165          do j=2,nbp_lat-1
166           ig0= 1+(j-2)*nbp_lon
167           do i=1,nbp_lon
168              area(i,j)=areafi_glo(ig0+i)
169           enddo
170           ! handle redundant point in longitude
171           area(nbp_lon+1,j)=area(1,j)
172          enddo
173         ENDIF
174
175         ! write "header" of file (longitudes, latitudes, geopotential, ...)
176         IF (klon_glo>1) THEN ! general 3D case
177           call iniwrite_specVI(nid,day_ini,area,nbp_lon+1,nbp_lat)
178         ELSE
179           call iniwrite_specVI(nid,day_ini,areafi_glo(1),1,1)
180         ENDIF
181         endif ! of if (is_master)
182
183         zitau = -1 ! initialize zitau
184      else
185         if (is_master) then
186           ! Open the NetCDF file
187           ierr = NF_OPEN(fichnom,NF_WRITE,nid)
188         endif
189      endif ! if (firstnom.eq.'1234567890')
190
191! Increment time index 'zitau' if it is the "firstcall" (at given time level)
192! to writediagfi
193!------------------------------------------------------------------------
194      if (nom.eq.firstnom) then
[3749]195          zitau = zitau + 1
[3184]196      end if
197
198!--------------------------------------------------------
199! Write the variables to output file if it's time to do so
200!--------------------------------------------------------
201
[3749]202      if ( MOD(zitau+1,isample) .eq.0.) then
[3184]203
204! Compute/write/extend 'Time' coordinate (date given in days)
205! (done every "first call" (at given time level) to writediagfi)
206! Note: date is incremented as 1 step ahead of physics time
207!       (like the 'histoire' outputs)
208!--------------------------------------------------------
209
210        if (nom.eq.firstnom) then
211
212        ! We have identified a "first call" (at given date)
213           ntime=ntime+1 ! increment # of stored time steps
214           ! compute corresponding date (in days and fractions thereof)
[3749]215           date= float (zitau +1)*(dtphys/daysec)
[3184]216
217           if (is_master) then
218             ! Get NetCDF ID of 'Time' variable
219             ierr= NF_INQ_VARID(nid,"Time",varid)
220
221             ! Write (append) the new date to the 'Time' array
222#ifdef NC_DOUBLE
223             ierr= NF_PUT_VARA_DOUBLE(nid,varid,ntime,1,date)
224#else
225             ierr= NF_PUT_VARA_REAL(nid,varid,ntime,1,date)
226#endif
227             if (ierr.ne.NF_NOERR) then
228              write(*,*) "***** PUT_VAR matter in writediagspec_nc"
229              write(*,*) "***** with time"
[3749]230              write(*,*) 'ierr=', ierr
[3184]231c             call abort
232             endif
233
234             write(6,*)'WRITEDIAGSPECVI: date= ', date
235           endif ! of if (is_master)
236        end if ! of if (nom.eq.firstnom)
237
238
[3749]239
[3184]240!Case of a 3D variable
241!---------------------
242        if (dimpx.eq.3) then
243
244!         A. Recast (copy) variable from physics grid to dynamics grid
245#ifdef CPP_PARA
246  ! gather field on a "global" (without redundant longitude) array
247          call Gather(px,dx3_glop)
248!$OMP MASTER
249          if (is_mpi_root) then
250            call Grid1Dto2D_glo(dx3_glop,dx3_glo)
251            ! copy dx3_glo() to dx3(:) and add redundant longitude
252            dx3(1:nbp_lon,:,:)=dx3_glo(1:nbp_lon,:,:)
253            dx3(nbp_lon+1,:,:)=dx3(1,:,:)
254          endif
255!$OMP END MASTER
256!$OMP BARRIER
257#else
258          IF (klon_glo>1) THEN ! General case
259           DO l=1,L_NSPECTV
260             DO i=1,nbp_lon+1
261                dx3(i,1,l)=px(1,l)
262                dx3(i,nbp_lat,l)=px(ngrid,l)
263             ENDDO
264             DO j=2,nbp_lat-1
265                ig0= 1+(j-2)*nbp_lon
266                DO i=1,nbp_lon
267                   dx3(i,j,l)=px(ig0+i,l)
268                ENDDO
269                dx3(nbp_lon+1,j,l)=dx3(1,j,l)
270             ENDDO
271           ENDDO
272          ELSE ! 1D model case
273           dx3_1d(1,1:L_NSPECTV)=px(1,1:L_NSPECTV)
274          ENDIF
275#endif
276
277!         B. Write (append) the variable to the NetCDF file
278          if (is_master) then
279
280! name of the variable
281           ierr= NF_INQ_VARID(nid,nom,varid)
282           if (ierr /= NF_NOERR) then
283! corresponding dimensions
284              ierr= NF_INQ_DIMID(nid,"longitude",id(1))
285              ierr= NF_INQ_DIMID(nid,"latitude",id(2))
286              ierr= NF_INQ_DIMID(nid,"VI_Wavenumber",id(3))
287              ierr= NF_INQ_DIMID(nid,"Time",id(4))
288
289! Create the variable if it doesn't exist yet
290
291              write (*,*) "=========================="
292              write (*,*) "DIAGSPEC: creating variable ",nom
293              call def_var(nid,nom,titre,unite,4,id,varid,ierr)
294
295           endif
296
297           corner(1)=1
298           corner(2)=1
299           corner(3)=1
300           corner(4)=ntime
301
302           IF (klon_glo==1) THEN
303             edges(1)=1
304           ELSE
305             edges(1)=nbp_lon+1
306           ENDIF
307           edges(2)=nbp_lat
308           edges(3)=L_NSPECTV
309           edges(4)=1
310#ifdef NC_DOUBLE
311           IF (klon_glo>1) THEN ! General case
312             ierr= NF_PUT_VARA_DOUBLE(nid,varid,corner,edges,dx3)
313           ELSE
314             ierr= NF_PUT_VARA_DOUBLE(nid,varid,corner,edges,dx3_1d)
315           ENDIF
316#else
317           IF (klon_glo>1) THEN ! General case
318             ierr= NF_PUT_VARA_REAL(nid,varid,corner,edges,dx3)
319           ELSE
320             ierr= NF_PUT_VARA_REAL(nid,varid,corner,edges,dx3_1d)
321           ENDIF
322#endif
323
324           if (ierr.ne.NF_NOERR) then
325              write(*,*) "***** PUT_VAR problem in writediagspec"
326              write(*,*) "***** with ",nom
327              write(*,*) 'ierr=', ierr
328             call abort
[3749]329           endif
[3184]330
331          endif ! of if (is_master)
332
333        endif ! of if (dimpx.eq.3)
334
335      endif ! of if ( MOD(zitau+1,irythme) .eq.0.)
336
337      ! Close the NetCDF file
338      if (is_master) then
339        ierr= NF_CLOSE(nid)
340      endif
341
342      end
Note: See TracBrowser for help on using the repository browser.