source: trunk/LMDZ.MARS/libf/phymars/physdem1.F @ 616

Last change on this file since 616 was 410, checked in by acolaitis, 14 years ago

Modified NCDF norm of our files from classic to 64bit offset to support variables indices of more than integer*4 maximum length. This is a priori retrocompatible with classic format. I have tested it, it works for all file outputs, newstart (on classic or 64-bit offset files). One can check the format of his .nc with ncdump -k file.nc
File size: 21.1 KB
RevLine 
[38]1      subroutine physdem1(filename,lonfi,latfi,nsoil,nq,
2     .                   phystep,day_ini,
3     .                   time,tsurf,tsoil,co2ice,emis,q2,qsurf,
4     .                   airefi,alb,ith,pzmea,pzstd,pzsig,pzgam,pzthe)
5
6      implicit none
7c-------------------------------------------------------------
8c
9c create physics (re-)start data file "restartfi.nc"
10c
11c
12c
13#include "dimensions.h"
14#include "paramet.h"
15c-----------------------------------------------------------------------
16#include "comvert.h"
17#include "comgeom2.h"
18#include "control.h"
19#include "comdissnew.h"
20#include "logic.h"
21#include "ener.h"
22#include "netcdf.inc"
23#include "dimphys.h"
24#include"advtrac.h"
25#include"callkeys.h"
26c
27      INTEGER nid,iq
28      INTEGER, parameter :: ivap=1
29      REAL, parameter :: qsolmax= 150.0
30      character (len=*) :: filename
31      character (len=7) :: str7
32
33      REAL day_ini
34      INTEGER nsoil,nq
35      integer ierr,idim1,idim2,idim3,idim4,idim5,nvarid
36
37c
38      REAL phystep,time
39      REAL latfi(ngridmx), lonfi(ngridmx)
40!      REAL champhys(ngridmx)
41      REAL tsurf(ngridmx)
42      INTEGER length
43      PARAMETER (length=100)
44      REAL tab_cntrl(length)
45
46c
47
48!      EXTERNAL defrun_new,iniconst,geopot,inigeom,massdair,pression
49!      EXTERNAL exner_hyb , SSUM
50c
51#include "serre.h"
52#include "clesph0.h"
53#include "fxyprim.h"
54#include "comgeomfi.h"
55#include "surfdat.h"
56#include "comsoil.h"
57#include "planete.h"
58#include "dimradmars.h"
59#include "yomaer.h"
60#include "comcstfi.h"
61
62      real co2ice(ngridmx),tsoil(ngridmx,nsoil),emis(ngridmx)
63      real q2(ngridmx, llm+1),qsurf(ngridmx,nq)
64      real airefi(ngridmx)
65      real alb(ngridmx),ith(ngridmx,nsoil)
66      real pzmea(ngridmx),pzstd(ngridmx)
67      real pzsig(ngridmx),pzgam(ngridmx),pzthe(ngridmx)
68      integer ig
69
70! flag which identifies if we are using old tracer names (qsurf01,...)
71      logical :: oldtracernames=.false.
72      integer :: count
73      character(len=30) :: txt ! to store some text
74! indexes of water vapour & water ice tracers (if any):
75      integer :: i_h2o_vap=0
76      integer :: i_h2o_ice=0
77c-----------------------------------------------------------------------
78
79      ! copy airefi(:) to area(:)
80      CALL SCOPY(ngridmx,airefi,1,area,1)
81      ! note: area() is defined in comgeomfi.h
82
83      DO ig=1,ngridmx
84         albedodat(ig)=alb(ig) ! note: albedodat() is defined in surfdat.h
85         zmea(ig)=pzmea(ig) ! note: zmea() is defined in surfdat.h
86         zstd(ig)=pzstd(ig) ! note: zstd() is defined in surfdat.h
87         zsig(ig)=pzsig(ig) ! note: zsig() is defined in surfdat.h
88         zgam(ig)=pzgam(ig) ! note: zgam() is defined in surfdat.h
89         zthe(ig)=pzthe(ig) ! note: zthe() is defined in surfdat.h
90      ENDDO
91
92      inertiedat(:,:)=ith(:,:) ! note inertiedat() is defined in comsoil.h
93c
94c  things to store in the physics start file:
95c
[410]96      ierr = NF_CREATE(adjustl(filename), 
97     &      IOR(NF_CLOBBER,NF_64BIT_OFFSET), nid)
[38]98      IF (ierr.NE.NF_NOERR) THEN
99        WRITE(6,*)'physdem1: Problem creating file ',adjustl(filename)
100        write(6,*) NF_STRERROR(ierr)
101        CALL ABORT
102      ENDIF
103c
104      ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 18,
105     .                       "Physics start file")
106c
107      ierr = NF_DEF_DIM (nid,"index",length,idim1)
108      if (ierr.ne.NF_NOERR) then
109        WRITE(6,*)'physdem1: Problem defining index dimension'
110        write(6,*) NF_STRERROR(ierr)
111        call abort
112      endif
113c
114      ierr = NF_DEF_DIM (nid,"physical_points",ngridmx,idim2)
115      if (ierr.ne.NF_NOERR) then
116        WRITE(6,*)'physdem1: Problem defining physical_points dimension'
117        write(6,*) NF_STRERROR(ierr)
118        call abort
119      endif
120c
121      ierr = NF_DEF_DIM (nid,"subsurface_layers",nsoil,idim3)
122      if (ierr.ne.NF_NOERR) then
123      WRITE(6,*)'physdem1: Problem defining subsurface_layers dimension'
124        write(6,*) NF_STRERROR(ierr)
125        call abort
126      endif
127c
128!      ierr = NF_DEF_DIM (nid,"nlayer+1",llm+1,idim4)
129      ierr = NF_DEF_DIM (nid,"nlayer_plus_1",llm+1,idim4)
130      if (ierr.ne.NF_NOERR) then
131        WRITE(6,*)'physdem1: Problem defining nlayer+1 dimension'
132        write(6,*) NF_STRERROR(ierr)
133        call abort
134      endif
135c
136      ierr = NF_DEF_DIM (nid,"number_of_advected_fields",nq,idim5)
137      if (ierr.ne.NF_NOERR) then
138        WRITE(6,*)'physdem1: Problem defining advected fields dimension'
139        WRITE(6,*)' nq = ',nq,' and ierr = ', ierr
140        write(6,*) NF_STRERROR(ierr)
141      endif
142
143      ierr = NF_ENDDEF(nid) ! exit NetCDF define mode
144
145c clear tab_cntrl(:) array
146      DO ierr = 1, length
147         tab_cntrl(ierr) = 0.0
148      ENDDO
149
150      write(*,*) "physdem1: ngridmx: ",ngridmx
151
152ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
153c Fill control array tab_cntrl(:) with paramleters for this run
154ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
155c Informations on the physics grid
156      tab_cntrl(1) = float(ngridmx)  ! number of nodes on physics grid
157      tab_cntrl(2) = float(nlayermx) ! number of atmospheric layers
158      tab_cntrl(3) = day_ini + int(time)         ! initial day
159      tab_cntrl(4) = time -int(time)            ! initiale time of day
160
161c Informations about Mars, used by dynamics and physics
162      tab_cntrl(5) = rad      ! radius of Mars (m) ~3397200
163      tab_cntrl(6) = omeg     ! rotation rate (rad.s-1)
164      tab_cntrl(7) = g        ! gravity (m.s-2) ~3.72
165      tab_cntrl(8) = mugaz    ! Molar mass of the atmosphere (g.mol-1) ~43.49
166      tab_cntrl(9) = rcp      !  = r/cp  ~0.256793 (=kappa dans dynamique)
167      tab_cntrl(10) = daysec  ! length of a sol (s)  ~88775
168
169      tab_cntrl(11) = phystep  ! time step in the physics
170      tab_cntrl(12) = 0.
171      tab_cntrl(13) = 0.
172
173c Informations about Mars, only for physics
174      tab_cntrl(14) = year_day  ! length of year (sols) ~668.6
175      tab_cntrl(15) = periheli  ! min. Sun-Mars distance (Mkm) ~206.66
176      tab_cntrl(16) = aphelie   ! max. SUn-Mars distance (Mkm) ~249.22
177      tab_cntrl(17) = peri_day  ! date of perihelion (sols since N. spring)
178      tab_cntrl(18) = obliquit  ! Obliquity of the planet (deg) ~23.98
179
180c Boundary layer and turbulence
[224]181      tab_cntrl(19) = z0_default   ! default surface roughness (m) ~0.01
[38]182      tab_cntrl(20) = lmixmin   ! mixing length ~100
183      tab_cntrl(21) = emin_turb ! minimal energy ~1.e-8
184
185c Optical properties of polar caps and ground emissivity
186      tab_cntrl(22) = albedice(1)  ! Albedo of northern cap ~0.5
187      tab_cntrl(23) = albedice(2)  ! Albedo of southern cap ~0.5
188      tab_cntrl(24) = emisice(1)   ! Emissivity of northern cap ~0.95
189      tab_cntrl(25) = emisice(2)   ! Emissivity of southern cap ~0.95
190      tab_cntrl(26) = emissiv      ! Emissivity of martian soil ~.95
191      tab_cntrl(31) = iceradius(1) ! mean scat radius of CO2 snow (north)
192      tab_cntrl(32) = iceradius(2) ! mean scat radius of CO2 snow (south)
193      tab_cntrl(33) = dtemisice(1) ! time scale for snow metamorphism (north)
194      tab_cntrl(34) = dtemisice(2) ! time scale for snow metamorphism (south)
195
196c dust aerosol properties
197      tab_cntrl(27) = tauvis      ! mean visible optical depth
198
199      tab_cntrl(28) = 0.
200      tab_cntrl(29) = 0.
201      tab_cntrl(30) = 0.
202
203! Soil properties:
204      tab_cntrl(35) = volcapa ! soil volumetric heat capacity
205     
206c
207!      write(*,*) "physdem1: tab_cntrl():",tab_cntrl
208     
209      ierr = NF_REDEF (nid) ! Enter NetCDF (re-)define mode
210      IF (ierr.NE.NF_NOERR) THEN
211         PRINT*, 'physdem1: Failed to swich to NetCDF define mode'
212         CALL abort
213      ENDIF
214      ! define variable
215#ifdef NC_DOUBLE
216      ierr = NF_DEF_VAR (nid, "controle", NF_DOUBLE, 1, idim1,nvarid)
217#else
218      ierr = NF_DEF_VAR (nid, "controle", NF_FLOAT, 1, idim1,nvarid)
219#endif
220      IF (ierr.NE.NF_NOERR) THEN
221         PRINT*, 'physdem1: Failed to define controle'
222         CALL abort
223      ENDIF
224      ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 18,
225     .                        "Control parameters")
226      IF (ierr.NE.NF_NOERR) THEN
227         PRINT*, 'physdem1: Failed to define controle title attribute'
228         CALL abort
229      ENDIF
230      ierr = NF_ENDDEF(nid) ! Leave NetCDF define mode
231      IF (ierr.NE.NF_NOERR) THEN
232         PRINT*, 'physdem1: Failed to swich out of NetCDF define mode'
233         CALL abort
234      ENDIF
235      ! write variable
236#ifdef NC_DOUBLE
237      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tab_cntrl)
238#else
239      ierr = NF_PUT_VAR_REAL (nid,nvarid,tab_cntrl)
240#endif
241      IF (ierr.NE.NF_NOERR) THEN
242         PRINT*, 'physdem1: Failed to write controle data'
243         CALL abort
244      ENDIF
245
246! write mid-layer depths mlayer() !known from comsoil.h
247
248      ierr = NF_REDEF (nid) ! Enter NetCDF (re-)define mode
249      ! define variable
250#ifdef NC_DOUBLE
251      ierr = NF_DEF_VAR (nid,"soildepth",NF_DOUBLE,1,idim3,nvarid)
252#else
253      ierr = NF_DEF_VAR (nid,"soildepth",NF_FLOAT,1,idim3,nvarid)
254#endif
255      ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 20,
256     .                        "Soil mid-layer depth")
257      ierr = NF_ENDDEF(nid) ! Leave NetCDF define mode
258      ! write variable
259#ifdef NC_DOUBLE
260      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,mlayer)
261#else
262      ierr = NF_PUT_VAR_REAL (nid,nvarid,mlayer)
263#endif
264
265c
266
267      ierr = NF_REDEF (nid)
268#ifdef NC_DOUBLE
269      ierr = NF_DEF_VAR (nid, "longitude", NF_DOUBLE, 1, idim2,nvarid)
270#else
271      ierr = NF_DEF_VAR (nid, "longitude", NF_FLOAT, 1, idim2,nvarid)
272#endif
273      ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 26,
274     .                        "Longitudes of physics grid")
275      ierr = NF_ENDDEF(nid)
276
277#ifdef NC_DOUBLE
278      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,lonfi)
279#else
280      ierr = NF_PUT_VAR_REAL (nid,nvarid,lonfi)
281#endif
282
283c
284
285      ierr = NF_REDEF (nid)
286#ifdef NC_DOUBLE
287      ierr = NF_DEF_VAR (nid, "latitude", NF_DOUBLE, 1, idim2,nvarid)
288#else
289      ierr = NF_DEF_VAR (nid, "latitude", NF_FLOAT, 1, idim2,nvarid)
290#endif
291      ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 25,
292     .                        "Latitudes of physics grid")
293      ierr = NF_ENDDEF(nid)
294#ifdef NC_DOUBLE
295      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,latfi)
296#else
297      ierr = NF_PUT_VAR_REAL (nid,nvarid,latfi)
298#endif
299
300c
301
302      ierr = NF_REDEF (nid)
303#ifdef NC_DOUBLE
304      ierr = NF_DEF_VAR (nid, "area", NF_DOUBLE, 1, idim2,nvarid)
305#else
306      ierr = NF_DEF_VAR (nid, "area", NF_FLOAT, 1, idim2,nvarid)
307#endif
308      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 9,
309     .                        "Mesh area")
310      ierr = NF_ENDDEF(nid)
311#ifdef NC_DOUBLE
312      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,area)
313#else
314      ierr = NF_PUT_VAR_REAL (nid,nvarid,area)
315#endif
316
317c
318
319      ierr = NF_REDEF (nid)
320#ifdef NC_DOUBLE
321      ierr = NF_DEF_VAR (nid, "phisfi", NF_DOUBLE, 1, idim2,nvarid)
322#else
323      ierr = NF_DEF_VAR (nid, "phisfi", NF_FLOAT, 1, idim2,nvarid)
324#endif
325      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 27,
326     .                        "Geopotential at the surface")
327      ierr = NF_ENDDEF(nid)
328#ifdef NC_DOUBLE
329      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,phisfi)
330#else
331      ierr = NF_PUT_VAR_REAL (nid,nvarid,phisfi)
332#endif
333
334c
335
336      ierr = NF_REDEF (nid)
337#ifdef NC_DOUBLE
338      ierr = NF_DEF_VAR (nid, "albedodat", NF_DOUBLE, 1, idim2,nvarid)
339#else
340      ierr = NF_DEF_VAR (nid, "albedodat", NF_FLOAT, 1, idim2,nvarid)
341#endif
342      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 21,
343     .                        "Albedo of bare ground")
344      ierr = NF_ENDDEF(nid)
345#ifdef NC_DOUBLE
346      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,albedodat)
347#else
348      ierr = NF_PUT_VAR_REAL (nid,nvarid,albedodat)
349#endif
350
351c
352c   some data for Francois Lott's programs
353c
354
355      ierr = NF_REDEF (nid)
356#ifdef NC_DOUBLE
357      ierr = NF_DEF_VAR (nid, "ZMEA", NF_DOUBLE, 1, idim2,nvarid)
358#else
359      ierr = NF_DEF_VAR (nid, "ZMEA", NF_FLOAT, 1, idim2,nvarid)
360#endif
361      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 19,
362     .                        "Relief: mean relief")
363      ierr = NF_ENDDEF(nid)
364#ifdef NC_DOUBLE
365      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zmea)
366#else
367      ierr = NF_PUT_VAR_REAL (nid,nvarid,zmea)
368#endif
369c
370      ierr = NF_REDEF (nid)
371#ifdef NC_DOUBLE
372      ierr = NF_DEF_VAR (nid, "ZSTD", NF_DOUBLE, 1, idim2,nvarid)
373#else
374      ierr = NF_DEF_VAR (nid, "ZSTD", NF_FLOAT, 1, idim2,nvarid)
375#endif
376      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 26,
377     .                        "Relief: standard deviation")
378      ierr = NF_ENDDEF(nid)
379#ifdef NC_DOUBLE
380      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zstd)
381#else
382      ierr = NF_PUT_VAR_REAL (nid,nvarid,zstd)
383#endif
384c
385      ierr = NF_REDEF (nid)
386#ifdef NC_DOUBLE
387      ierr = NF_DEF_VAR (nid, "ZSIG", NF_DOUBLE, 1, idim2,nvarid)
388#else
389      ierr = NF_DEF_VAR (nid, "ZSIG", NF_FLOAT, 1, idim2,nvarid)
390#endif
391      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 23,
392     .                        "Relief: sigma parameter")
393      ierr = NF_ENDDEF(nid)
394#ifdef NC_DOUBLE
395      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zsig)
396#else
397      ierr = NF_PUT_VAR_REAL (nid,nvarid,zsig)
398#endif
399c
400      ierr = NF_REDEF (nid)
401#ifdef NC_DOUBLE
402      ierr = NF_DEF_VAR (nid, "ZGAM", NF_DOUBLE, 1, idim2,nvarid)
403#else
404      ierr = NF_DEF_VAR (nid, "ZGAM", NF_FLOAT, 1, idim2,nvarid)
405#endif
406      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 23,
407     .                        "Relief: gamma parameter")
408      ierr = NF_ENDDEF(nid)
409#ifdef NC_DOUBLE
410      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zgam)
411#else
412      ierr = NF_PUT_VAR_REAL (nid,nvarid,zgam)
413#endif
414c
415      ierr = NF_REDEF (nid)
416#ifdef NC_DOUBLE
417      ierr = NF_DEF_VAR (nid, "ZTHE", NF_DOUBLE, 1, idim2,nvarid)
418#else
419      ierr = NF_DEF_VAR (nid, "ZTHE", NF_FLOAT, 1, idim2,nvarid)
420#endif
421      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 23,
422     .                        "Relief: theta parameter")
423      ierr = NF_ENDDEF(nid)
424#ifdef NC_DOUBLE
425      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,zthe)
426#else
427      ierr = NF_PUT_VAR_REAL (nid,nvarid,zthe)
428#endif
429
430c Write the physical fields
431
432! CO2 Ice Cover
433
434      ierr = NF_REDEF (nid)
435#ifdef NC_DOUBLE
436      ierr = NF_DEF_VAR (nid, "co2ice", NF_DOUBLE, 1, idim2,nvarid)
437#else
438      ierr = NF_DEF_VAR (nid, "co2ice", NF_FLOAT, 1, idim2,nvarid)
439#endif
440      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 13,
441     .                        "CO2 ice cover")
442      ierr = NF_ENDDEF(nid)
443#ifdef NC_DOUBLE
444      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,co2ice)
445#else
446      ierr = NF_PUT_VAR_REAL (nid,nvarid,co2ice)
447#endif
448
449! Soil Thermal inertia
450
451      ierr = NF_REDEF (nid)
452#ifdef NC_DOUBLE
453      ierr = NF_DEF_VAR (nid,"inertiedat",NF_DOUBLE,
454     &                   2,(/idim2,idim3/),nvarid)
455#else
456      ierr = NF_DEF_VAR (nid,"inertiedat",NF_FLOAT,
457     &                   2,(/idim2,idim3/),nvarid)
458#endif
459      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 20,
460     .                        "Soil thermal inertia")
461      ierr = NF_ENDDEF(nid)
462#ifdef NC_DOUBLE
463      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,inertiedat)
464#else
465      ierr = NF_PUT_VAR_REAL (nid,nvarid,inertiedat)
466#endif
467
468!  Surface temperature
469
470      ierr = NF_REDEF (nid)
471#ifdef NC_DOUBLE
472      ierr = NF_DEF_VAR (nid, "tsurf", NF_DOUBLE, 1, idim2,nvarid)
473#else
474      ierr = NF_DEF_VAR (nid, "tsurf", NF_FLOAT, 1, idim2,nvarid)
475#endif
476      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 19,
477     .                        "Surface temperature")
478      ierr = NF_ENDDEF(nid)
479#ifdef NC_DOUBLE
480      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tsurf)
481#else
482      ierr = NF_PUT_VAR_REAL (nid,nvarid,tsurf)
483#endif
484
485! Soil temperature
486
487      ierr = NF_REDEF (nid)
488#ifdef NC_DOUBLE
489      ierr = NF_DEF_VAR (nid,"tsoil",NF_DOUBLE,2,(/idim2,idim3/),nvarid)
490#else
491!      ierr = NF_DEF_VAR (nid, "tsoil", NF_FLOAT, 2, idim2,nvarid)
492      ierr = NF_DEF_VAR (nid,"tsoil",NF_FLOAT,2,(/idim2,idim3/),nvarid)
493#endif
494      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 16,
495     .                        "Soil temperature")
496      ierr = NF_ENDDEF(nid)
497#ifdef NC_DOUBLE
498      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tsoil)
499#else
500      ierr = NF_PUT_VAR_REAL (nid,nvarid,tsoil)
501#endif
502
503c emissivity
504
505      ierr = NF_REDEF (nid)
506#ifdef NC_DOUBLE
507      ierr = NF_DEF_VAR (nid, "emis", NF_DOUBLE, 1, idim2,nvarid)
508#else
509      ierr = NF_DEF_VAR (nid, "emis", NF_FLOAT, 1, idim2,nvarid)
510#endif
511      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 18,
512     .                        "Surface emissivity")
513      ierr = NF_ENDDEF(nid)
514#ifdef NC_DOUBLE
515      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,emis)
516#else
517      ierr = NF_PUT_VAR_REAL (nid,nvarid,emis)
518#endif
519
[224]520! surface roughness length (z0 is a common in surfdat.h)
521
522      ierr = NF_REDEF (nid)
523#ifdef NC_DOUBLE
524      ierr = NF_DEF_VAR (nid, "z0", NF_DOUBLE, 1, idim2,nvarid)
525#else
526      ierr = NF_DEF_VAR (nid, "z0", NF_FLOAT, 1, idim2,nvarid)
527#endif
528      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 24,
529     .                        "Surface roughness length")
530      ierr = NF_ENDDEF(nid)
531#ifdef NC_DOUBLE
532      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,z0)
533#else
534      ierr = NF_PUT_VAR_REAL (nid,nvarid,z0)
535#endif
536
537
[38]538c planetary boundary layer
539
540      ierr = NF_REDEF (nid)
541#ifdef NC_DOUBLE
542      ierr = NF_DEF_VAR (nid,"q2",NF_DOUBLE,2,(/idim2,idim4/),nvarid)
543#else
544      ierr = NF_DEF_VAR (nid,"q2",NF_FLOAT, 2,(/idim2,idim4/),nvarid)
545#endif
546      ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 17,
547     .                        "pbl wind variance")
548      ierr = NF_ENDDEF(nid)
549#ifdef NC_DOUBLE
550      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,q2)
551#else
552      ierr = NF_PUT_VAR_REAL (nid,nvarid,q2)
553#endif
554      IF (ierr.NE.NF_NOERR) THEN
555        PRINT*, 'physdem1: Failed to write q2'
556        CALL abort
557      ENDIF
558
559c tracers
560
561! Preliminary stuff: check if tracers follow old naming convention (qsurf01,
562!                    qsurf02, ...)
563      count=0
564      do iq=1,nqmx
565        txt= " "
566        write(txt,'(a1,i2.2)')'q',iq
567        if (txt.ne.tnom(iq)) then ! use tracer names stored in dynamics
568          ! did not find old tracer name
569          exit ! might as well stop here
570        else
571          ! found old tracer name
572          count=count+1
573        endif
574      enddo
575      if (count.eq.nqmx) then
576        write(*,*) "physdem1:tracers seem to follow old naming ",
577     &             "convention (qsurf01,qsurf02,...)"
578        write(*,*) "   => will work for now ... "
579        write(*,*) "      but you should run newstart to rename them"
580        oldtracernames=.true.
581        ! Moreover, if computing water cycle with ice, move surface ice
582        ! back to qsurf(nqmx)
583        IF (water) THEN
[164]584          !"loop" to avoid potential out-of-bounds on arrays
[38]585          write(*,*)'physdem1: moving surface water ice to index ',nqmx
[164]586          do iq=nqmx,nqmx
587          qsurf(1:ngridmx,iq)=qsurf(1:ngridmx,iq-1)
588          qsurf(1:ngridmx,iq-1)=0
589          enddo
[38]590        ENDIF
591      endif ! of if (count.eq.nqmx)
592
593      IF(nq.GE.1) THEN
594! preliminary stuff: look for water vapour & water ice tracers (if any)
595        if (.not.oldtracernames) then
596         do iq=1,nq
597           if (tnom(iq).eq."h2o_vap") then
598             i_h2o_vap=iq
599           endif
600           if (tnom(iq).eq."h2o_ice") then
601             i_h2o_ice=iq
602           endif
603         enddo ! of iq=1,nq
604         ! handle special case of only water vapour tracer (no ice)
605         if ((i_h2o_vap.ne.0).and.(i_h2o_ice.eq.0)) then
606          ! then the index of (surface) ice is i_h2o_vap
607          i_h2o_ice=i_h2o_vap
608         endif
609        endif ! of if (.not.oldtracernames)
610
611         DO iq=1,nq
612           IF (oldtracernames) THEN
613             txt=" "
614             write(txt,'(a5,i2.2)')'qsurf',iq
615           ELSE
616             txt=tnom(iq)
617             ! Exception: there is no water vapour surface tracer
618             if (txt.eq."h2o_vap") then
619               write(*,*)"physdem1: skipping water vapour tracer"
620               if (i_h2o_ice.eq.i_h2o_vap) then
621               ! then there is no "water ice" tracer; but still
622               ! there is some water ice on the surface
623                 write(*,*)"          writting water ice instead"
624                 txt="h2o_ice"
625               else
626               ! there is a "water ice" tracer which has been / will be
627               ! delt with in due time
628                 cycle
629               endif ! of if (igcm_h2o_ice.eq.igcm_h2o_vap)
630             endif ! of if (txt.eq."h2o_vap")
631           ENDIF ! of IF (oldtracernames)
632
633           ierr=NF_REDEF(nid)
634           IF (ierr.NE.NF_NOERR) THEN
635             PRINT*, 'physdem1: Failed to swich to NetCDF define mode'
636             CALL abort
637           ENDIF
638#ifdef NC_DOUBLE
639           ierr=NF_DEF_VAR(nid,txt,NF_DOUBLE,1,idim2,nvarid)
640#else
641           ierr=NF_DEF_VAR(nid,txt,NF_FLOAT,1,idim2,nvarid)
642#endif
643           IF (ierr.NE.NF_NOERR) THEN
644             PRINT*, 'physdem1: Failed to define ',trim(txt)
645             CALL abort
646           ENDIF
647           ierr=NF_PUT_ATT_TEXT (nid, nvarid, "title", 17,
648     &                        "tracer on surface")
649           IF (ierr.NE.NF_NOERR) THEN
650             PRINT*, 'physdem1: Failed to define ',trim(txt),
651     &               ' title attribute'
652             CALL abort
653           ENDIF
654           ierr=NF_ENDDEF(nid)
655           IF (ierr.NE.NF_NOERR) THEN
656             PRINT*, 'physdem1: Failed to swich out of define mode'
657             CALL abort
658           ENDIF
659           
660#ifdef NC_DOUBLE
661            ierr=NF_PUT_VAR_DOUBLE (nid,nvarid,qsurf(1,iq))
662#else
663            ierr=NF_PUT_VAR_REAL (nid,nvarid,qsurf(1,iq))
664#endif
665           IF (ierr.NE.NF_NOERR) THEN
666             PRINT*, 'physdem1: Failed to write ',trim(txt)
667             CALL abort
668           ENDIF
669         ENDDO ! of DO iq=1,nq
670      ENDIF ! of IF(nq.GE.1)
671
672c close file
673      ierr = NF_CLOSE(nid)
674
675      RETURN
676
677      END
Note: See TracBrowser for help on using the repository browser.