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

Last change on this file since 3992 was 3928, checked in by jmauxion, 3 months ago

Generic PCM:
Adding a slow_diagfi flag to the run.def/rcm1d.def file for 1D models only. When False, the netcdf
file is opened/closed once, thus saving significant computing time. When true,
the opening frequency is at output frequency (recommended in debug mode). Also
fixing a redundant loop on tracers when writing outputs in physiq_mod.
JM

File size: 11.6 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: diagfi_output_rate, dtphys, daysec
51      use time_phylmdz_mod, only: day_ini
52      use callkeys_mod, only: iradia
53
54      implicit none
55
56      include "netcdf.inc"
57
58! Arguments on input:
59      integer,intent(in) :: ngrid
60      character (len=*),intent(in) :: nom,titre,unite
61      integer,intent(in) :: dimpx
62      real,intent(in) :: px(ngrid,L_NSPECTI)
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
70      integer isample
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 OLR output
91      real dx3(nbp_lon+1,nbp_lat,L_NSPECTI) ! to store the data set
92      real dx3_1d(1,L_NSPECTI) ! 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_NSPECTI)
97      real dx3_glo(nbp_lon,nbp_lat,L_NSPECTI) ! 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
102      if (grid_type == unstructured) then
103        return
104      endif
105
106!***************************************************************
107!Sortie des variables au rythme voulu
108
109      isample = diagfi_output_rate*iradia ! output rate
110!EM+JL if the spectra need to be output more frequently, need to define a
111! diagspec_output_rate...
112
113
114!***************************************************************
115
116! Initialisation of 'firstnom' and create/open the "diagfi.nc" NetCDF file
117! ------------------------------------------------------------------------
118! (Au tout premier appel de la subroutine durant le run.)
119
120      fichnom="diagspecIR.nc"
121
122      if (firstnom.eq.'1234567890') then ! .true. for the very first call
123      !  to this subroutine; now set 'firstnom'
124         firstnom = nom
125         ! just to be sure, check that firstnom is large enough to hold nom
126         if (len_trim(firstnom).lt.len_trim(nom)) then
127           write(*,*) "writediagspecIR: Error !!!"
128           write(*,*) "   firstnom string not long enough!!"
129           write(*,*) "   increase its size to at least ",len_trim(nom)
130           call abort_physic("writediagspecIR",
131     &          "firstnom is not large enough to hold nom" ,1)
132         endif
133
134#ifdef CPP_PARA
135          ! Gather cell_area() mesh area on physics grid
136          call Gather(cell_area,areafi_glo)
137#else
138         areafi_glo(:)=cell_area(:)
139#endif
140         ! Create the NetCDF file
141         if (is_master) then
142         ierr = NF_CREATE(fichnom, NF_CLOBBER, nid)
143         ! Define the 'Time' dimension
144         ierr = nf_def_dim(nid,"Time",NF_UNLIMITED,idim)
145         ! Define the 'Time' variable
146#ifdef NC_DOUBLE
147         ierr = NF_DEF_VAR (nid, "Time", NF_DOUBLE, 1, idim,varid)
148#else
149         ierr = NF_DEF_VAR (nid, "Time", NF_FLOAT, 1, idim,varid)
150#endif
151         ! Add a long_name attribute
152         ierr = NF_PUT_ATT_TEXT (nid, varid, "long_name",
153     .          4,"Time")
154         ! Add a units attribute
155         ierr = NF_PUT_ATT_TEXT(nid, varid,'units',29,
156     .          "days since 0000-00-0 00:00:00")
157         ! Switch out of NetCDF Define mode
158         ierr = NF_ENDDEF(nid)
159
160         ! Build area()
161         IF (klon_glo>1) THEN
162          do i=1,nbp_lon+1 ! poles
163           ! divide at the poles by nbp_lon
164           area(i,1)=areafi_glo(1)/nbp_lon
165           area(i,nbp_lat)=areafi_glo(klon_glo)/nbp_lon
166          enddo
167          do j=2,nbp_lat-1
168           ig0= 1+(j-2)*nbp_lon
169           do i=1,nbp_lon
170              area(i,j)=areafi_glo(ig0+i)
171           enddo
172           ! handle redundant point in longitude
173           area(nbp_lon+1,j)=area(1,j)
174          enddo
175         ENDIF
176
177         ! write "header" of file (longitudes, latitudes, area, ...)
178         IF (klon_glo>1) THEN ! general 3D case
179           call iniwrite_specIR(nid,day_ini,area,nbp_lon+1,nbp_lat)
180         ELSE
181           call iniwrite_specIR(nid,day_ini,areafi_glo(1),1,1)
182         ENDIF
183           ! Close the NetCDF file
184           ierr= NF_CLOSE(nid)
185         endif ! of if (is_master)
186
187         zitau = -1 ! initialize zitau
188      endif ! if (firstnom.eq.'1234567890')
189
190! Increment time index 'zitau' if it is the "firstcall" (at given time level)
191! to writediagfi
192!------------------------------------------------------------------------
193      if (nom.eq.firstnom) then
194          zitau = zitau + 1
195      end if
196
197!--------------------------------------------------------
198! Write the variables to output file if it's time to do so
199!--------------------------------------------------------
200
201      if ( MOD(zitau+1,isample) .eq.0.) then
202
203         if (is_master) then
204           ! Open the NetCDF file
205           ierr = NF_OPEN(fichnom,NF_WRITE,nid)
206         endif
207
208! Compute/write/extend 'Time' coordinate (date given in days)
209! (done every "first call" (at given time level) to writediagfi)
210! Note: date is incremented as 1 step ahead of physics time
211!       (like the 'histoire' outputs)
212!--------------------------------------------------------
213
214        if (nom.eq.firstnom) then
215
216        ! We have identified a "first call" (at given date)
217           ntime=ntime+1 ! increment # of stored time steps
218           ! compute corresponding date (in days and fractions thereof)
219           date=(zitau +1)*(dtphys/daysec)
220
221           if (is_master) then
222             ! Get NetCDF ID of 'Time' variable
223             ierr= NF_INQ_VARID(nid,"Time",varid)
224
225             ! Write (append) the new date to the 'Time' array
226#ifdef NC_DOUBLE
227             ierr= NF_PUT_VARA_DOUBLE(nid,varid,ntime,1,date)
228#else
229             ierr= NF_PUT_VARA_REAL(nid,varid,ntime,1,date)
230#endif
231             if (ierr.ne.NF_NOERR) then
232              write(*,*) "***** PUT_VAR matter in writediagspec_nc"
233              write(*,*) "***** with time"
234              write(*,*) 'ierr=', ierr   
235c             call abort
236             endif
237
238             write(6,*)'WRITEDIAGSPECIR: date= ', date
239           endif ! of if (is_master)
240        end if ! of if (nom.eq.firstnom)
241
242
243 
244!Case of a 3D variable
245!---------------------
246        if (dimpx.eq.3) then
247
248!         A. Recast (copy) variable from physics grid to dynamics grid
249#ifdef CPP_PARA
250  ! gather field on a "global" (without redundant longitude) array
251          call Gather(px,dx3_glop)
252!$OMP MASTER
253          if (is_mpi_root) then
254            call Grid1Dto2D_glo(dx3_glop,dx3_glo)
255            ! copy dx3_glo() to dx3(:) and add redundant longitude
256            dx3(1:nbp_lon,:,:)=dx3_glo(1:nbp_lon,:,:)
257            dx3(nbp_lon+1,:,:)=dx3(1,:,:)
258          endif
259!$OMP END MASTER
260!$OMP BARRIER
261#else
262          IF (klon_glo>1) THEN ! General case
263           DO l=1,L_NSPECTI
264             DO i=1,nbp_lon+1
265                dx3(i,1,l)=px(1,l)
266                dx3(i,nbp_lat,l)=px(ngrid,l)
267             ENDDO
268             DO j=2,nbp_lat-1
269                ig0= 1+(j-2)*nbp_lon
270                DO i=1,nbp_lon
271                   dx3(i,j,l)=px(ig0+i,l)
272                ENDDO
273                dx3(nbp_lon+1,j,l)=dx3(1,j,l)
274             ENDDO
275           ENDDO
276          ELSE ! 1D model case
277            dx3_1d(1,1:L_NSPECTI)=px(1,1:L_NSPECTI)
278          ENDIF
279#endif
280
281!         B. Write (append) the variable to the NetCDF file
282          if (is_master) then
283
284! name of the variable
285           ierr= NF_INQ_VARID(nid,nom,varid)
286           if (ierr /= NF_NOERR) then
287! corresponding dimensions
288              ierr= NF_INQ_DIMID(nid,"longitude",id(1))
289              ierr= NF_INQ_DIMID(nid,"latitude",id(2))
290              ierr= NF_INQ_DIMID(nid,"IR_Wavenumber",id(3))
291              ierr= NF_INQ_DIMID(nid,"Time",id(4))
292
293! Create the variable if it doesn't exist yet
294
295              write (*,*) "=========================="
296              write (*,*) "DIAGSPECIR: creating variable ",nom
297              call def_var(nid,nom,titre,unite,4,id,varid,ierr)
298
299           endif
300
301           corner(1)=1
302           corner(2)=1
303           corner(3)=1
304           corner(4)=ntime
305
306           IF (klon_glo==1) THEN
307             edges(1)=1
308           ELSE
309             edges(1)=nbp_lon+1
310           ENDIF
311           edges(2)=nbp_lat
312           edges(3)=L_NSPECTI
313           edges(4)=1
314#ifdef NC_DOUBLE
315           IF (klon_glo>1) THEN ! General case
316             ierr= NF_PUT_VARA_DOUBLE(nid,varid,corner,edges,dx3)
317           ELSE
318             ierr= NF_PUT_VARA_DOUBLE(nid,varid,corner,edges,dx3_1d)
319           ENDIF
320#else
321           IF (klon_glo>1) THEN ! General case
322             ierr= NF_PUT_VARA_REAL(nid,varid,corner,edges,dx3)
323           ELSE
324             ierr= NF_PUT_VARA_REAL(nid,varid,corner,edges,dx3_1d)
325           ENDIF
326#endif
327
328           if (ierr.ne.NF_NOERR) then
329              write(*,*) "***** PUT_VAR problem in writediagspec"
330              write(*,*) "***** with ",nom
331              write(*,*) 'ierr=', ierr
332             call abort_physic("diagspecIR",
333     &            "PUT_VAR problem in writediagspec", 1)
334           endif
335
336          endif ! of if (is_master)
337
338        endif ! of if (dimpx.eq.3)
339
340        ! Close the NetCDF file
341        if (is_master) then
342          ierr= NF_CLOSE(nid)
343        endif
344
345      endif ! of if ( MOD(zitau+1,isample) .eq.0.)
346
347      end
Note: See TracBrowser for help on using the repository browser.