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

Last change on this file since 3787 was 3725, checked in by emillour, 3 months ago

Generic PCM:
Change the way the rate of outputs for diagfi.nc files is specified:
IMPORTANT: Specifying "ecritphy" no longer possible and will trigger an error.
Use "diagfi_output_rate" to specify output rate (in physics steps) instead.
This should makes things (hopefully) clearer for users and also better
enforces a cleaner and clearer separation between dynamics and physics.
EM

File size: 11.3 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           stop
131         endif
132
133#ifdef CPP_PARA
134          ! Gather cell_area() mesh area on physics grid
135          call Gather(cell_area,areafi_glo)
136#else
137         areafi_glo(:)=cell_area(:)
138#endif
139         ! Create the NetCDF file
140         if (is_master) then
141         ierr = NF_CREATE(fichnom, NF_CLOBBER, nid)
142         ! Define the 'Time' dimension
143         ierr = nf_def_dim(nid,"Time",NF_UNLIMITED,idim)
144         ! Define the 'Time' variable
145#ifdef NC_DOUBLE
146         ierr = NF_DEF_VAR (nid, "Time", NF_DOUBLE, 1, idim,varid)
147#else
148         ierr = NF_DEF_VAR (nid, "Time", NF_FLOAT, 1, idim,varid)
149#endif
150         ! Add a long_name attribute
151         ierr = NF_PUT_ATT_TEXT (nid, varid, "long_name",
152     .          4,"Time")
153         ! Add a units attribute
154         ierr = NF_PUT_ATT_TEXT(nid, varid,'units',29,
155     .          "days since 0000-00-0 00:00:00")
156         ! Switch out of NetCDF Define mode
157         ierr = NF_ENDDEF(nid)
158
159         ! Build area()
160         IF (klon_glo>1) THEN
161          do i=1,nbp_lon+1 ! poles
162           ! divide at the poles by nbp_lon
163           area(i,1)=areafi_glo(1)/nbp_lon
164           area(i,nbp_lat)=areafi_glo(klon_glo)/nbp_lon
165          enddo
166          do j=2,nbp_lat-1
167           ig0= 1+(j-2)*nbp_lon
168           do i=1,nbp_lon
169              area(i,j)=areafi_glo(ig0+i)
170           enddo
171           ! handle redundant point in longitude
172           area(nbp_lon+1,j)=area(1,j)
173          enddo
174         ENDIF
175
176         ! write "header" of file (longitudes, latitudes, area, ...)
177         IF (klon_glo>1) THEN ! general 3D case
178           call iniwrite_specIR(nid,day_ini,area,nbp_lon+1,nbp_lat)
179         ELSE
180           call iniwrite_specIR(nid,day_ini,areafi_glo(1),1,1)
181         ENDIF
182         endif ! of if (is_master)
183
184         zitau = -1 ! initialize zitau
185      else
186         if (is_master) then
187           ! Open the NetCDF file
188           ierr = NF_OPEN(fichnom,NF_WRITE,nid)
189         endif
190      endif ! if (firstnom.eq.'1234567890')
191
192! Increment time index 'zitau' if it is the "firstcall" (at given time level)
193! to writediagfi
194!------------------------------------------------------------------------
195      if (nom.eq.firstnom) then
196          zitau = zitau + 1
197      end if
198
199!--------------------------------------------------------
200! Write the variables to output file if it's time to do so
201!--------------------------------------------------------
202
203      if ( MOD(zitau+1,isample) .eq.0.) then
204
205! Compute/write/extend 'Time' coordinate (date given in days)
206! (done every "first call" (at given time level) to writediagfi)
207! Note: date is incremented as 1 step ahead of physics time
208!       (like the 'histoire' outputs)
209!--------------------------------------------------------
210
211        if (nom.eq.firstnom) then
212
213        ! We have identified a "first call" (at given date)
214           ntime=ntime+1 ! increment # of stored time steps
215           ! compute corresponding date (in days and fractions thereof)
216           date=(zitau +1)*(dtphys/daysec)
217
218           if (is_master) then
219             ! Get NetCDF ID of 'Time' variable
220             ierr= NF_INQ_VARID(nid,"Time",varid)
221
222             ! Write (append) the new date to the 'Time' array
223#ifdef NC_DOUBLE
224             ierr= NF_PUT_VARA_DOUBLE(nid,varid,ntime,1,date)
225#else
226             ierr= NF_PUT_VARA_REAL(nid,varid,ntime,1,date)
227#endif
228             if (ierr.ne.NF_NOERR) then
229              write(*,*) "***** PUT_VAR matter in writediagspec_nc"
230              write(*,*) "***** with time"
231              write(*,*) 'ierr=', ierr   
232c             call abort
233             endif
234
235             write(6,*)'WRITEDIAGSPECIR: date= ', date
236           endif ! of if (is_master)
237        end if ! of if (nom.eq.firstnom)
238
239
240 
241!Case of a 3D variable
242!---------------------
243        if (dimpx.eq.3) then
244
245!         A. Recast (copy) variable from physics grid to dynamics grid
246#ifdef CPP_PARA
247  ! gather field on a "global" (without redundant longitude) array
248          call Gather(px,dx3_glop)
249!$OMP MASTER
250          if (is_mpi_root) then
251            call Grid1Dto2D_glo(dx3_glop,dx3_glo)
252            ! copy dx3_glo() to dx3(:) and add redundant longitude
253            dx3(1:nbp_lon,:,:)=dx3_glo(1:nbp_lon,:,:)
254            dx3(nbp_lon+1,:,:)=dx3(1,:,:)
255          endif
256!$OMP END MASTER
257!$OMP BARRIER
258#else
259          IF (klon_glo>1) THEN ! General case
260           DO l=1,L_NSPECTI
261             DO i=1,nbp_lon+1
262                dx3(i,1,l)=px(1,l)
263                dx3(i,nbp_lat,l)=px(ngrid,l)
264             ENDDO
265             DO j=2,nbp_lat-1
266                ig0= 1+(j-2)*nbp_lon
267                DO i=1,nbp_lon
268                   dx3(i,j,l)=px(ig0+i,l)
269                ENDDO
270                dx3(nbp_lon+1,j,l)=dx3(1,j,l)
271             ENDDO
272           ENDDO
273          ELSE ! 1D model case
274            dx3_1d(1,1:L_NSPECTI)=px(1,1:L_NSPECTI)
275          ENDIF
276#endif
277
278!         B. Write (append) the variable to the NetCDF file
279          if (is_master) then
280
281! name of the variable
282           ierr= NF_INQ_VARID(nid,nom,varid)
283           if (ierr /= NF_NOERR) then
284! corresponding dimensions
285              ierr= NF_INQ_DIMID(nid,"longitude",id(1))
286              ierr= NF_INQ_DIMID(nid,"latitude",id(2))
287              ierr= NF_INQ_DIMID(nid,"IR_Wavenumber",id(3))
288              ierr= NF_INQ_DIMID(nid,"Time",id(4))
289
290! Create the variable if it doesn't exist yet
291
292              write (*,*) "=========================="
293              write (*,*) "DIAGSPECIR: creating variable ",nom
294              call def_var(nid,nom,titre,unite,4,id,varid,ierr)
295
296           endif
297
298           corner(1)=1
299           corner(2)=1
300           corner(3)=1
301           corner(4)=ntime
302
303           IF (klon_glo==1) THEN
304             edges(1)=1
305           ELSE
306             edges(1)=nbp_lon+1
307           ENDIF
308           edges(2)=nbp_lat
309           edges(3)=L_NSPECTI
310           edges(4)=1
311#ifdef NC_DOUBLE
312           IF (klon_glo>1) THEN ! General case
313             ierr= NF_PUT_VARA_DOUBLE(nid,varid,corner,edges,dx3)
314           ELSE
315             ierr= NF_PUT_VARA_DOUBLE(nid,varid,corner,edges,dx3_1d)
316           ENDIF
317#else
318           IF (klon_glo>1) THEN ! General case
319             ierr= NF_PUT_VARA_REAL(nid,varid,corner,edges,dx3)
320           ELSE
321             ierr= NF_PUT_VARA_REAL(nid,varid,corner,edges,dx3_1d)
322           ENDIF
323#endif
324
325           if (ierr.ne.NF_NOERR) then
326              write(*,*) "***** PUT_VAR problem in writediagspec"
327              write(*,*) "***** with ",nom
328              write(*,*) 'ierr=', ierr
329             call abort
330           endif
331
332          endif ! of if (is_master)
333
334        endif ! of if (dimpx.eq.3)
335
336      endif ! of if ( MOD(zitau+1,isample) .eq.0.)
337
338      ! Close the NetCDF file
339      if (is_master) then
340        ierr= NF_CLOSE(nid)
341      endif
342
343      end
Note: See TracBrowser for help on using the repository browser.