source: trunk/LMDZ.GENERIC/libf/phystd/phyetat0.F @ 704

Last change on this file since 704 was 699, checked in by emillour, 13 years ago

Generic GCM:

  • Corrected the polar mesh surface area which was wrong in the physics (changes in phyetat0.F, calfis.F and newstart.F)
  • Some cleanup in newstart.F (removed some obsolete "Mars" options: mons_ice,..) and also added option "q=profile" to initialize a tracer with a profile read from file "profile_tracername"

EM

File size: 19.8 KB
Line 
1      SUBROUTINE phyetat0 (fichnom,tab0,Lmodif,nsoil,nq,
2     .           day_ini,time,
3     .           tsurf,tsoil,emis,q2,qsurf,cloudfrac,totcloudfrac,hice)
4      implicit none
5c======================================================================
6c Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
7c  Adaptation à Mars : Yann Wanherdrick
8c Objet: Lecture de l etat initial pour la physique
9c======================================================================
10#include "netcdf.inc"
11#include "dimensions.h"
12#include "dimphys.h"
13!#include "comgeomfi.h"
14#include "surfdat.h"
15#include "planete.h"
16#include "comcstfi.h"
17#include"advtrac.h"
18
19c======================================================================
20      INTEGER nbsrf !Mars nbsrf a 1 au lieu de 4
21      PARAMETER (nbsrf=1) ! nombre de sous-fractions pour une maille
22!======================================================================
23!  Arguments:
24!  ---------
25!  inputs:
26      character*(*) fichnom ! "startfi.nc" file
27      integer tab0
28      integer Lmodif
29      integer nsoil ! # of soil layers
30      integer nq
31      integer day_ini
32      real time
33
34!  outputs:
35      real tsurf(ngridmx,nbsrf) ! surface temperature
36      real tsoil(ngridmx,nsoil,nbsrf) ! soil temperature
37      real emis(ngridmx) ! surface emissivity
38      real q2(ngridmx, llm+1) !
39      real qsurf(ngridmx,nq) ! tracers on surface
40!      real co2ice(ngridmx) ! co2 ice cover
41      real cloudfrac(ngridmx,nlayermx)
42      real hice(ngridmx), totcloudfrac(ngridmx)
43
44
45
46!======================================================================
47!  Local variables:
48
49!      INTEGER radpas
50!      REAL co2_ppm
51!      REAL solaire
52
53      real xmin,xmax ! to display min and max of a field
54c
55      INTEGER ig,iq,lmax
56      INTEGER nid, nvarid
57      INTEGER ierr, i, nsrf
58!      integer isoil
59!      INTEGER length
60!      PARAMETER (length=100)
61      CHARACTER*7 str7
62      CHARACTER*2 str2
63      CHARACTER*1 yes
64c
65      REAL p_rad,p_omeg,p_g,p_cpp,p_mugaz,p_daysec
66      INTEGER nqold
67
68! flag which identifies if 'startfi.nc' file is using old names (qsurf01,...)
69      logical :: oldtracernames=.false.
70      integer :: count
71      character(len=30) :: txt ! to store some text
72
73 
74c
75c Ouvrir le fichier contenant l etat initial:
76c
77
78      ierr = NF_OPEN (fichnom, NF_NOWRITE,nid)
79      IF (ierr.NE.NF_NOERR) THEN
80        write(6,*)' Pb d''ouverture du fichier '//fichnom
81        CALL ABORT
82      ENDIF
83
84! Preliminary stuff: check if tracers follow old naming convention (qsurf01,
85!                    qsurf02, ...)
86      count=0
87      do iq=1,nqmx
88        txt= " "
89        write(txt,'(a5,i2.2)')'qsurf',iq
90        ierr=NF_INQ_VARID(nid,txt,nvarid)
91        if (ierr.ne.NF_NOERR) then
92          ! did not find old tracer name
93          exit ! might as well stop here
94        else
95          ! found old tracer name
96          count=count+1
97        endif
98      enddo
99      if (count.eq.nqmx) then
100        write(*,*) "phyetat0:tracers seem to follow old naming ",
101     &             "convention (qsurf01,qsurf02,...)"
102        write(*,*) "   => will work for now ... "
103        write(*,*) "      but you should run newstart to rename them"
104        oldtracernames=.true.
105      endif
106
107c modifications possibles des variables de tab_cntrl
108      PRINT*
109      write(*,*) 'TABFI de phyeta0',Lmodif,tab0
110      call tabfi (nid,Lmodif,tab0,day_ini,lmax,p_rad,
111     .              p_omeg,p_g,p_cpp,p_mugaz,p_daysec,time)
112c
113c Read latitudes (coordinates): No need, these are provided by the dynamics
114c
115!      ierr=nf90_inq_varid(nid,"latitude",nvarid)
116!      IF (ierr.NE.nf90_noerr) THEN
117!         PRINT*, 'phyetat0: Le champ <latitude> est absent'
118!         write(*,*)trim(nf90_strerror(ierr))
119!         CALL abort
120!      ENDIF
121!      ierr=nf90_get_var(nid,nvarid,lati)
122!      IF (ierr.NE.nf90_noerr) THEN
123!         PRINT*, 'phyetat0: Lecture echouee pour <latitude>'
124!         write(*,*)trim(nf90_strerror(ierr))
125!         CALL abort
126!      ENDIF
127c
128c read longitudes (coordinates): No need, these are provided by the dynamics
129c
130!      ierr=nf90_inq_varid(nid,"longitude",nvarid)
131!      IF (ierr.NE.nf90_noerr) THEN
132!         PRINT*, 'phyetat0: Le champ <longitude> est absent'
133!         write(*,*)trim(nf90_strerror(ierr))
134!         CALL abort
135!      ENDIF
136!      ierr=nf90_get_var(nid,nvarid,long)
137!      IF (ierr.NE.nf90_noerr) THEN
138!         PRINT*, 'phyetat0: Lecture echouee pour <longitude>'
139!         write(*,*)trim(nf90_strerror(ierr))
140!         CALL abort
141!      ENDIF
142c
143c Read areas of meshes: No need, these are provided by the dynamics
144c
145!      ierr=nf90_inq_varid(nid,"area",nvarid)
146!      IF (ierr.NE.nf90_noerr) THEN
147!         PRINT*, 'phyetat0: Le champ <area> est absent'
148!         write(*,*)trim(nf90_strerror(ierr))
149!         CALL abort
150!      ENDIF
151!      ierr=nf90_get_var(nid,nvarid,area)
152!      IF (ierr.NE.nf90_noerr) THEN
153!         PRINT*, 'phyetat0: Lecture echouee pour <area>'
154!         write(*,*)trim(nf90_strerror(ierr))
155!         CALL abort
156!      ENDIF
157!      xmin = 1.0E+20
158!      xmax = -1.0E+20
159!      xmin = MINVAL(area)
160!      xmax = MAXVAL(area)
161!      PRINT*,'Aires des mailles <area>:', xmin, xmax
162c
163c Lecture du geopotentiel au sol:
164c
165      ierr = NF_INQ_VARID (nid, "phisfi", nvarid)
166      IF (ierr.NE.NF_NOERR) THEN
167         PRINT*, 'phyetat0: Le champ <phisfi> est absent'
168         CALL abort
169      ENDIF
170#ifdef NC_DOUBLE
171      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, phisfi)
172#else
173      ierr = NF_GET_VAR_REAL(nid, nvarid, phisfi)
174#endif
175      IF (ierr.NE.NF_NOERR) THEN
176         PRINT*, 'phyetat0: Lecture echouee pour <phisfi>'
177         CALL abort
178      ENDIF
179      xmin = 1.0E+20
180      xmax = -1.0E+20
181      xmin = MINVAL(phisfi)
182      xmax = MAXVAL(phisfi)
183      PRINT*,'Geopotentiel au sol <phisfi>:', xmin, xmax
184c
185c Lecture de l''albedo du sol nu:
186c
187      ierr = NF_INQ_VARID (nid, "albedodat", nvarid)
188      IF (ierr.NE.NF_NOERR) THEN
189         PRINT*, 'phyetat0: Le champ <albedodat> est absent'
190         CALL abort
191      ENDIF
192#ifdef NC_DOUBLE
193      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, albedodat)
194#else
195      ierr = NF_GET_VAR_REAL(nid, nvarid, albedodat)
196#endif
197      IF (ierr.NE.NF_NOERR) THEN
198         PRINT*, 'phyetat0: Lecture echouee pour <albedodat>'
199         CALL abort
200      ENDIF
201      xmin = 1.0E+20
202      xmax = -1.0E+20
203      xmin = MINVAL(albedodat)
204      xmax = MAXVAL(albedodat)
205      PRINT*,'Albedo du sol nu <albedodat>:', xmin, xmax
206c
207c Lecture de l''inertie thermique du sol:
208c
209!      ierr = NF_INQ_VARID (nid, "inertiedat", nvarid)
210!      IF (ierr.NE.NF_NOERR) THEN
211!         PRINT*, 'phyetat0: Le champ <inertiedat> est absent'
212!         CALL abort
213!      ENDIF
214!#ifdef NC_DOUBLE
215!      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, inertiedat)
216!#else
217!      ierr = NF_GET_VAR_REAL(nid, nvarid, inertiedat)
218!#endif
219!      IF (ierr.NE.NF_NOERR) THEN
220!         PRINT*, 'phyetat0: Lecture echouee pour <inertiedat>'
221!         CALL abort
222!      ENDIF
223!      xmin = 1.0E+20
224!      xmax = -1.0E+20
225!      xmin = MINVAL(inertiedat)
226!      xmax = MAXVAL(inertiedat)
227!      PRINT*,'Inertie thermique du sol <inertiedat>:', xmin, xmax
228c
229c ZMEA
230c
231      ierr = NF_INQ_VARID (nid, "ZMEA", nvarid)
232      IF (ierr.NE.NF_NOERR) THEN
233         PRINT*, 'phyetat0: Le champ <ZMEA> est absent'
234         CALL abort
235      ENDIF
236#ifdef NC_DOUBLE
237      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, zmea)
238#else
239      ierr = NF_GET_VAR_REAL(nid, nvarid, zmea)
240#endif
241      IF (ierr.NE.NF_NOERR) THEN
242         PRINT*, 'phyetat0: Lecture echouee pour <ZMEA>'
243         CALL abort
244      ENDIF
245      xmin = 1.0E+20
246      xmax = -1.0E+20
247      DO i = 1, ngridmx
248         xmin = MIN(zmea(i),xmin)
249         xmax = MAX(zmea(i),xmax)
250      ENDDO
251      PRINT*,'<zmea>:', xmin, xmax
252c
253c ZSTD
254c
255      ierr = NF_INQ_VARID (nid, "ZSTD", nvarid)
256      IF (ierr.NE.NF_NOERR) THEN
257         PRINT*, 'phyetat0: Le champ <ZSTD> est absent'
258         CALL abort
259      ENDIF
260#ifdef NC_DOUBLE
261      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, zstd)
262#else
263      ierr = NF_GET_VAR_REAL(nid, nvarid, zstd)
264#endif
265      IF (ierr.NE.NF_NOERR) THEN
266         PRINT*, 'phyetat0: Lecture echouee pour <ZSTD>'
267         CALL abort
268      ENDIF
269      xmin = 1.0E+20
270      xmax = -1.0E+20
271      DO i = 1, ngridmx
272         xmin = MIN(zstd(i),xmin)
273         xmax = MAX(zstd(i),xmax)
274      ENDDO
275      PRINT*,'<zstd>:', xmin, xmax
276c
277c ZSIG
278c
279      ierr = NF_INQ_VARID (nid, "ZSIG", nvarid)
280      IF (ierr.NE.NF_NOERR) THEN
281         PRINT*, 'phyetat0: Le champ <ZSIG> est absent'
282         CALL abort
283      ENDIF
284#ifdef NC_DOUBLE
285      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, zsig)
286#else
287      ierr = NF_GET_VAR_REAL(nid, nvarid, zsig)
288#endif
289      IF (ierr.NE.NF_NOERR) THEN
290         PRINT*, 'phyetat0: Lecture echouee pour <ZSIG>'
291         CALL abort
292      ENDIF
293      xmin = 1.0E+20
294      xmax = -1.0E+20
295      DO i = 1, ngridmx
296         xmin = MIN(zsig(i),xmin)
297         xmax = MAX(zsig(i),xmax)
298      ENDDO
299      PRINT*,'<zsig>:', xmin, xmax
300c
301c ZGAM
302c
303      ierr = NF_INQ_VARID (nid, "ZGAM", nvarid)
304      IF (ierr.NE.NF_NOERR) THEN
305         PRINT*, 'phyetat0: Le champ <ZGAM> est absent'
306         CALL abort
307      ENDIF
308#ifdef NC_DOUBLE
309      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, zgam)
310#else
311      ierr = NF_GET_VAR_REAL(nid, nvarid, zgam)
312#endif
313      IF (ierr.NE.NF_NOERR) THEN
314         PRINT*, 'phyetat0: Lecture echouee pour <ZGAM>'
315         CALL abort
316      ENDIF
317      xmin = 1.0E+20
318      xmax = -1.0E+20
319      DO i = 1, ngridmx
320         xmin = MIN(zgam(i),xmin)
321         xmax = MAX(zgam(i),xmax)
322      ENDDO
323      PRINT*,'<zgam>:', xmin, xmax
324c
325c ZTHE
326c
327      ierr = NF_INQ_VARID (nid, "ZTHE", nvarid)
328      IF (ierr.NE.NF_NOERR) THEN
329         PRINT*, 'phyetat0: Le champ <ZTHE> est absent'
330         CALL abort
331      ENDIF
332#ifdef NC_DOUBLE
333      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, zthe)
334#else
335      ierr = NF_GET_VAR_REAL(nid, nvarid, zthe)
336#endif
337      IF (ierr.NE.NF_NOERR) THEN
338         PRINT*, 'phyetat0: Lecture echouee pour <ZTHE>'
339         CALL abort
340      ENDIF
341      xmin = 1.0E+20
342      xmax = -1.0E+20
343      DO i = 1, ngridmx
344         xmin = MIN(zthe(i),xmin)
345         xmax = MAX(zthe(i),xmax)
346      ENDDO
347      PRINT*,'<zthe>:', xmin, xmax
348c
349c CO2 ice cover
350c
351! Ehouarn: from now on, there is no "co2ice" standalone field; it is supposed
352! to be with stored in qsurf(i_co2_ice)
353!      ierr = NF_INQ_VARID (nid, "co2ice", nvarid)
354!      IF (ierr.NE.NF_NOERR) THEN
355!         PRINT*, 'phyetat0: Le champ <co2ice> est absent'
356!         CALL abort
357!      ENDIF
358!#ifdef NC_DOUBLE
359!      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, co2ice)
360!#else
361!      ierr = NF_GET_VAR_REAL(nid, nvarid, co2ice)
362!#endif
363!      IF (ierr.NE.NF_NOERR) THEN
364!         PRINT*, 'phyetat0: Lecture echouee pour <co2ice>'
365!         CALL abort
366!      ENDIF
367!      xmin = 1.0E+20
368!      xmax = -1.0E+20
369!      xmin = MINVAL(co2ice)
370!      xmax = MAXVAL(co2ice)
371!      PRINT*,'CO2 ice cover <co2ice>:', xmin, xmax
372c
373c Lecture des temperatures du sol:
374c
375      ierr = NF_INQ_VARID (nid, "tsurf", nvarid)
376      IF (ierr.NE.NF_NOERR) THEN
377         PRINT*, 'phyetat0: Le champ <tsurf> est absent'
378         PRINT*, '          Mais je vais essayer de lire TS**'
379         IF (nbsrf.GT.99) THEN
380            PRINT*, "Trop de sous-mailles"
381            CALL abort
382         ENDIF
383         DO nsrf = 1, nbsrf
384           WRITE(str2,'(i2.2)') nsrf
385           ierr = NF_INQ_VARID (nid, "TS"//str2, nvarid)
386           IF (ierr.NE.NF_NOERR) THEN
387           PRINT*, "phyetat0: Le champ <TS"//str2//"> est absent"
388              CALL abort
389           ENDIF
390#ifdef NC_DOUBLE
391           ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tsurf(1,nsrf))
392#else
393           ierr = NF_GET_VAR_REAL(nid, nvarid, tsurf(1,nsrf))
394#endif
395           IF (ierr.NE.NF_NOERR) THEN
396             PRINT*, "phyetat0: Lecture echouee pour <TS"//str2//">"
397             CALL abort
398           ENDIF
399           xmin = 1.0E+20
400           xmax = -1.0E+20
401           xmin = MINVAL(tsurf)
402           xmax = MAXVAL(tsurf)
403           PRINT*,'Temperature du sol TS**:', nsrf, xmin, xmax
404         ENDDO
405      ELSE
406         PRINT*, 'phyetat0: Le champ <tsurf> est present'
407         PRINT*, '          J ignore donc les autres temperatures TS**'
408#ifdef NC_DOUBLE
409         ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tsurf(1,1))
410#else
411         ierr = NF_GET_VAR_REAL(nid, nvarid, tsurf(1,1))
412#endif
413         IF (ierr.NE.NF_NOERR) THEN
414            PRINT*, "phyetat0: Lecture echouee pour <TSURF>"
415            CALL abort
416         ENDIF
417         xmin = 1.0E+20
418         xmax = -1.0E+20
419         xmin = MINVAL(tsurf)
420         xmax = MAXVAL(tsurf)
421         PRINT*,'Temperature du sol <tsurf>', xmin, xmax
422         IF (nbsrf >= 2) THEN
423            DO nsrf = 2, nbsrf
424               DO i = 1, ngridmx
425                  tsurf(i,nsrf) = tsurf(i,1)
426               ENDDO
427            ENDDO
428         ENDIF
429      ENDIF
430c
431c Lecture des temperatures du sol profond:
432c
433!      IF (nsoil.GT.99 .OR. nbsrf.GT.99) THEN
434!         PRINT*, "Trop de couches ou sous-mailles"
435!         CALL abort
436!      ENDIF
437!      DO nsrf = 1, nbsrf
438!         DO isoil=1, nsoil
439!            WRITE(str7,'(i2.2,"srf",i2.2)') isoil, nsrf
440!            ierr = NF_INQ_VARID (nid, 'tsoil', nvarid)
441!            IF (ierr.NE.NF_NOERR) THEN
442!               PRINT*, "phyetat0: Le champ <tsoil> est absent"
443!               PRINT*, "          Il prend donc la valeur de surface"
444!               DO i=1, ngridmx
445!                  tsoil(i,isoil,nsrf)=tsurf(i,nsrf)
446!               ENDDO
447!            ELSE
448!#ifdef NC_DOUBLE
449!              ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tsoil(1,1,nsrf))
450!#else
451!              ierr = NF_GET_VAR_REAL(nid, nvarid, tsoil(1,1,nsrf))
452!#endif
453!               IF (ierr.NE.NF_NOERR) THEN
454!                  PRINT*, "Lecture echouee pour <tsoil>"
455!                  CALL abort
456!               ENDIF
457!            ENDIF
458!         ENDDO
459!      ENDDO
460!      xmin = 1.0E+20
461!      xmax = -1.0E+20
462!      xmin = MINVAL(tsoil)
463!      xmax = MAXVAL(tsoil)
464!      PRINT*,'Temperatures du sol profond <tsoil>', xmin, xmax
465c
466c Surface emissivity
467c
468      ierr = NF_INQ_VARID (nid, "emis", nvarid)
469      IF (ierr.NE.NF_NOERR) THEN
470         PRINT*, 'phyetat0: Le champ <emis> est absent'
471         CALL abort
472      ENDIF
473#ifdef NC_DOUBLE
474      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, emis)
475#else
476      ierr = NF_GET_VAR_REAL(nid, nvarid, emis)
477#endif
478      IF (ierr.NE.NF_NOERR) THEN
479         PRINT*, 'phyetat0: Lecture echouee pour <emis>'
480         CALL abort
481      ENDIF
482      xmin = 1.0E+20
483      xmax = -1.0E+20
484      xmin = MINVAL(emis)
485      xmax = MAXVAL(emis)
486      PRINT*,'Surface emissivity <emis>:', xmin, xmax
487
488c
489c Cloud fraction (added by BC 2010)
490c
491      ierr = NF_INQ_VARID (nid, "cloudfrac", nvarid)
492      IF (ierr.NE.NF_NOERR) THEN
493         PRINT*, 'phyetat0: Le champ <cloudfrac> est absent'
494      cloudfrac(:,:)=0.5
495!         CALL abort
496      ENDIF
497#ifdef NC_DOUBLE
498      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, cloudfrac)
499#else
500      ierr = NF_GET_VAR_REAL(nid, nvarid, cloudfrac)
501#endif
502      IF (ierr.NE.NF_NOERR) THEN
503         PRINT*, 'phyetat0: Lecture echouee pour <cloudfrac>'
504         CALL abort
505      ENDIF
506      xmin = 1.0E+20
507      xmax = -1.0E+20
508      xmin = MINVAL(cloudfrac)
509      xmax = MAXVAL(cloudfrac)
510      PRINT*,'Cloud fraction <cloudfrac>:', xmin, xmax
511
512
513c
514c Total cloud fraction (added by BC 2010)
515c
516      ierr = NF_INQ_VARID (nid, "totcloudfrac", nvarid)
517!      ierr = NF_INQ_VARID (nid, "totalfrac", nvarid)
518      IF (ierr.NE.NF_NOERR) THEN
519         PRINT*, 'phyetat0: Le champ <totcloudfrac> est absent'
520      totcloudfrac(:)=0.5
521!         CALL abort
522      ENDIF
523#ifdef NC_DOUBLE
524      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, totcloudfrac)
525#else
526      ierr = NF_GET_VAR_REAL(nid, nvarid, totcloudfrac)
527#endif
528      IF (ierr.NE.NF_NOERR) THEN
529         PRINT*, 'phyetat0: Lecture echouee pour <totcloudfrac>'
530         !CALL abort
531      ENDIF
532      xmin = 1.0E+20
533      xmax = -1.0E+20
534      xmin = MINVAL(totcloudfrac)
535      xmax = MAXVAL(totcloudfrac)
536      PRINT*,'Cloud fraction <totcloudfrac>:', xmin, xmax
537
538
539
540c
541c Height of oceanic ice (added by BC 2010)
542c
543      ierr = NF_INQ_VARID (nid, "hice", nvarid)
544      IF (ierr.NE.NF_NOERR) THEN
545         PRINT*, 'phyetat0: Le champ <hice> est absent'
546         CALL abort
547      ENDIF
548#ifdef NC_DOUBLE
549      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, hice)
550#else
551      ierr = NF_GET_VAR_REAL(nid, nvarid, hice)
552#endif
553      IF (ierr.NE.NF_NOERR) THEN
554         PRINT*, 'phyetat0: Lecture echouee pour <hice>'
555         CALL abort
556      ENDIF
557      xmin = 1.0E+20
558      xmax = -1.0E+20
559      xmin = MINVAL(hice)
560      xmax = MAXVAL(hice)
561      PRINT*,'Height of oceanic ice <hice>:', xmin, xmax
562
563
564c
565c pbl wind variance
566c
567      ierr = NF_INQ_VARID (nid, "q2", nvarid)
568      IF (ierr.NE.NF_NOERR) THEN
569         PRINT*, 'phyetat0: Le champ <q2> est absent'
570         CALL abort
571      ENDIF
572#ifdef NC_DOUBLE
573      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, q2)
574#else
575      ierr = NF_GET_VAR_REAL(nid, nvarid, q2)
576#endif
577      IF (ierr.NE.NF_NOERR) THEN
578         PRINT*, 'phyetat0: Lecture echouee pour <q2>'
579         CALL abort
580      ENDIF
581      xmin = 1.0E+20
582      xmax = -1.0E+20
583      xmin = MINVAL(q2)
584      xmax = MAXVAL(q2)
585      PRINT*,'pbl wind variance <q2>:', xmin, xmax
586c
587c tracer on surface
588c
589
590      IF(nq.GE.1) THEN
591         nqold=nq
592         DO iq=1,nq
593!            str7(1:5)='qsurf'
594!            WRITE(str7(6:7),'(i2.2)') iq
595!            ierr = NF_INQ_VARID (nid,str7,nvarid)
596           IF (oldtracernames) THEN
597             txt=" "
598             write(txt,'(a5,i2.2)')'qsurf',iq
599           ELSE
600             txt=tnom(iq)
601!             if (txt.eq."h2o_vap") then
602               ! There is no surface tracer for h2o_vap;
603               ! "h2o_ice" should be loaded instead
604!               txt="h2o_ice"
605!               write(*,*) 'phyetat0: loading surface tracer',
606!     &                     ' h2o_ice instead of h2o_vap'
607!             endif
608           ENDIF ! of IF (oldtracernames) THEN
609           ierr=NF_INQ_VARID(nid,txt,nvarid)
610           IF (ierr.NE.NF_NOERR) THEN
611             write(*,*) 'PHYETAT0: WARNING : surface tracer',trim(txt),
612     &                  ' not found in file'
613             write(*,*) trim(txt), ' set to 0'
614             do ig=1,ngridmx
615               qsurf(ig,iq)=0.
616             end do
617             nqold=min(iq-1,nqold)
618           ELSE
619#ifdef NC_DOUBLE
620             ierr = NF_GET_VAR_DOUBLE(nid, nvarid,qsurf(1,iq))
621#else
622             ierr = NF_GET_VAR_REAL(nid, nvarid,qsurf(1,iq))
623#endif
624             IF (ierr.NE.NF_NOERR) THEN
625               PRINT*, 'phyetat0: Lecture echouee pour <',trim(txt),'>'
626               CALL abort
627             ENDIF
628           ENDIF
629           xmin = 1.0E+20
630           xmax = -1.0E+20
631           xmin = MINVAL(qsurf(1:ngridmx,iq))
632           xmax = MAXVAL(qsurf(1:ngridmx,iq))
633           PRINT*,'tracer on surface <',trim(txt),'>:',xmin,xmax
634         ENDDO
635         if ((nqold.lt.nq).and.(nqold.ge.1)) then
636c        case when new tracer are added in addition to old ones
637             write(*,*)'qsurf 1 to ', nqold,'were already present'
638             write(*,*)'qsurf ', nqold+1,' to ', nqmx,'are new'
639             write(*,*)' and initialized to zero'
640             qsurf(:,nqold+1:nqmx)=0.0
641!            yes=' '
642!            do while ((yes.ne.'y').and.(yes.ne.'n'))
643!             write(*,*) 'Would you like to reindex qsurf # 1 ->',nqold
644!             write(*,*) 'to #',nqmx-nqold+1,'->', nqmx,'   (y or n) ?'
645!             read(*,fmt='(a)') yes
646!            end do
647!            if (yes.eq.'y') then
648!              write(*,*) 'OK, let s reindex qsurf'
649!                 do ig=1,ngridmx
650!                    do iq=nqmx,nqmx-nqold+1,-1
651!                       qsurf(ig,iq)=qsurf(ig,iq-nqmx+nqold)
652!                    end do
653!                    do iq=nqmx-nqold,1,-1
654!                       qsurf(ig,iq)= 0.
655!                    end do
656!                 end do
657!            end if
658         end if ! of if ((nqold.lt.nq).and.(nqold.ge.1))
659      ENDIF ! of IF(nq.GE.1)
660
661! Call to soil_settings, in order to read soil temperatures,
662! as well as thermal inertia and volumetric heat capacity
663
664      call soil_settings(nid,ngridmx,nsoil,tsurf,tsoil)
665c
666c Fermer le fichier:
667c
668      ierr = NF_CLOSE(nid)
669c
670      RETURN
671      END
Note: See TracBrowser for help on using the repository browser.