source: trunk/LMDZ.GENERIC/libf/phystd/writediagspecIR.F @ 3058

Last change on this file since 3058 was 2732, checked in by aslmd, 2 years ago

took into account unstructured grid in writediagspec*.F

File size: 11.4 KB
Line 
1      subroutine writediagspecIR(ngrid,nom,titre,unite,dimpx,px)
2
3!  Ecriture de variables diagnostiques au choix dans la physique
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)
8!  Dans la version 2000, la periode d'ecriture est celle de
9!  "ecritphy " regle dans le fichier de controle de run :  run.def
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
17! date donnee sont legerement differentes que dans le fichier histoire car
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
22!
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)
28!     
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!=================================================================
42 
43! Addition by RW (2010) to allow OLR to be saved in .nc format
44      use radinc_h, only : L_NSPECTI
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
50      use time_phylmdz_mod, only: ecritphy, iphysiq, day_step, day_ini
51      use callkeys_mod, only: iradia
52
53      implicit none
54
55      include "netcdf.inc"
56
57! Arguments on input:
58      integer,intent(in) :: ngrid
59      character (len=*),intent(in) :: nom,titre,unite
60      integer,intent(in) :: dimpx
61      real,intent(in) :: px(ngrid,L_NSPECTI)
62
63! Local variables:
64
65!      real dx3(iip1,jjp1,llm) ! to store a 3D data set
66!      real dx2(iip1,jjp1)     ! to store a 2D (surface) data set
67!      real dx0
68
69      integer irythme
70      integer ierr
71      integer iq
72      integer i,j,l,zmax , ig0
73
74      integer,save :: zitau=0
75      character(len=20),save :: firstnom='1234567890'
76      real,save :: date
77!$OMP THREADPRIVATE(firstnom,zitau,date)
78
79! Ajouts
80      integer, save :: ntime=0
81!$OMP THREADPRIVATE(ntime)
82      integer :: idim,varid
83      integer :: nid
84      character (len =50):: fichnom
85      integer, dimension(4) :: id
86      integer, dimension(4) :: edges,corner
87
88      real area((nbp_lon+1),nbp_lat)
89! added by RDW for OLR output
90      real dx3(nbp_lon+1,nbp_lat,L_NSPECTI) ! to store the data set
91      real dx3_1d(1,L_NSPECTI) ! to store the data with 1D model
92
93#ifdef CPP_PARA
94! Added to work in parallel mode
95      real dx3_glop(klon_glo,L_NSPECTI)
96      real dx3_glo(nbp_lon,nbp_lat,L_NSPECTI) ! to store a global 3D data set
97      real areafi_glo(klon_glo) ! mesh area on global physics grid
98#else
99      real areafi_glo(ngrid) ! mesh area on global physics grid
100#endif
101      if (grid_type == unstructured) then
102        return
103      endif
104
105!***************************************************************
106!Sortie des variables au rythme voulu
107
108      irythme = ecritphy*iradia ! sortie au rythme de ecritphy*iradia
109!EM+JL if the spetra need to be output more frequently, need to define a ecritSpec...
110!     irythme = iphysiq  ! sortie a tous les pas physique
111
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="diagspecIR.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(*,*) "writediagspecIR: 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, area, ...)
176         IF (klon_glo>1) THEN ! general 3D case
177           call iniwrite_specIR(nid,day_ini,area,nbp_lon+1,nbp_lat)
178         ELSE
179           call iniwrite_specIR(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
195          zitau = zitau + iphysiq
196      end if
197
198!--------------------------------------------------------
199! Write the variables to output file if it's time to do so
200!--------------------------------------------------------
201
202      if ( MOD(zitau+1,irythme) .eq.0.) then
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)
215           date= float (zitau +1)/float (day_step)
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"
230              write(*,*) 'ierr=', ierr   
231c             call abort
232             endif
233
234             write(6,*)'WRITEDIAGSPECIR: date= ', date
235           endif ! of if (is_master)
236        end if ! of if (nom.eq.firstnom)
237
238
239 
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_NSPECTI
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_NSPECTI)=px(1,1:L_NSPECTI)
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,"IR_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 (*,*) "DIAGSPECIR: 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_NSPECTI
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
329           endif
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.