source: trunk/LMDZ.GENERIC/libf/dyn3d/start2archive.F @ 588

Last change on this file since 588 was 253, checked in by emillour, 13 years ago

Generic GCM

  • Massive update to version 0.7

EM+RW

File size: 14.9 KB
Line 
1c=======================================================================
2      PROGRAM start2archive
3c=======================================================================
4c
5c
6c   Date:    01/1997
7c   ----
8c
9c
10c   Objet:   Passage des  fichiers netcdf d'etat initial "start" et
11c   -----    "startfi" a un fichier netcdf unique "start_archive"
12c
13c  "start_archive" est une banque d'etats initiaux:
14c  On peut stocker plusieurs etats initiaux dans un meme fichier "start_archive"
15c    (Veiller dans ce cas avoir un day_ini different pour chacun des start)
16c
17c
18c
19c=======================================================================
20
21      implicit none
22
23#include "dimensions.h"
24#include "paramet.h"
25#include "comconst.h"
26#include "comdissip.h"
27#include "comvert.h"
28#include "comgeom.h"
29#include "logic.h"
30#include "temps.h"
31#include "control.h"
32#include "ener.h"
33#include "description.h"
34
35#include "dimphys.h"
36#include "comsoil.h"
37#include"advtrac.h"
38#include "netcdf.inc"
39
40c-----------------------------------------------------------------------
41c   Declarations
42c-----------------------------------------------------------------------
43
44c variables dynamiques du GCM
45c -----------------------------
46      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm) ! vents covariants
47      REAL teta(ip1jmp1,llm)                    ! temperature potentielle
48      REAL q(ip1jmp1,llm,nqmx)               ! champs advectes
49      REAL pks(ip1jmp1)                      ! exner (f pour filtre)
50      REAL pk(ip1jmp1,llm)
51      REAL pkf(ip1jmp1,llm)
52      REAL beta(iip1,jjp1,llm)
53      REAL phis(ip1jmp1)                     ! geopotentiel au sol
54      REAL masse(ip1jmp1,llm)                ! masse de l'atmosphere
55      REAL ps(ip1jmp1)                       ! pression au sol
56      REAL p3d(iip1, jjp1, llm+1)            ! pression aux interfaces
57     
58c Variable Physiques (grille physique)
59c ------------------------------------
60      REAL tsurf(ngridmx)       ! Surface temperature
61      REAL tsoil(ngridmx,nsoilmx) ! Soil temperature
62      REAL co2ice(ngridmx)      ! CO2 ice layer
63      REAL q2(ngridmx,nlayermx+1),qsurf(ngridmx,nqmx)
64      REAL emis(ngridmx)
65      INTEGER start,length
66      PARAMETER (length = 100)
67      REAL tab_cntrl_fi(length) ! tableau des parametres de startfi
68      INTEGER*4 day_ini_fi
69
70!     added by FF for cloud fraction setup
71      REAL hice(ngridmx)
72      REAL cloudfrac(ngridmx,nlayermx),totalcloudfrac(ngridmx)
73
74
75c Variable naturelle / grille scalaire
76c ------------------------------------
77      REAL T(ip1jmp1,llm),us(ip1jmp1,llm),vs(ip1jmp1,llm)
78      REAL tsurfS(ip1jmp1)
79      REAL tsoilS(ip1jmp1,nsoilmx)
80      REAL ithS(ip1jmp1,nsoilmx) ! Soil Thermal Inertia
81      REAL co2iceS(ip1jmp1)
82      REAL q2S(ip1jmp1,llm+1),qsurfS(ip1jmp1,nqmx)
83      REAL emisS(ip1jmp1)
84
85!     added by FF for cloud fraction setup
86      REAL hiceS(ip1jmp1)
87      REAL cloudfracS(ip1jmp1,llm),totalcloudfracS(ip1jmp1)
88
89c Variables intermediaires : vent naturel, mais pas coord scalaire
90c----------------------------------------------------------------
91      REAL vn(ip1jm,llm),un(ip1jmp1,llm)
92
93c Autres  variables
94c -----------------
95      LOGICAL startdrs
96      INTEGER Lmodif
97
98      REAL ptotal, co2icetotal
99      REAL timedyn,timefi !fraction du jour dans start, startfi
100      REAL date
101
102      CHARACTER*2 str2
103      CHARACTER*80 fichier
104      data  fichier /'startfi'/
105
106      INTEGER ij, l,i,j,isoil,iq
107      character*80      fichnom
108      integer :: ierr,ntime
109      integer :: nq,numvanle
110      character(len=30) :: txt ! to store some text
111
112c Netcdf
113c-------
114      integer varid,dimid,timelen
115      INTEGER nid,nid1
116
117c-----------------------------------------------------------------------
118c   Initialisations
119c-----------------------------------------------------------------------
120
121      grireg   = .TRUE.
122
123c=======================================================================
124c Lecture des donnees
125c=======================================================================
126! Load tracer names:
127      call iniadvtrac(nq,numvanle)
128
129      fichnom = 'start.nc'
130      CALL dynetat0(fichnom,nqmx,vcov,ucov,teta,q,masse,
131     .       ps,phis,timedyn)
132
133
134      fichnom = 'startfi.nc'
135      Lmodif=0
136
137      CALL phyetat0 (fichnom,0,Lmodif,nsoilmx,nqmx,day_ini_fi,timefi,
138     .      tsurf,tsoil,emis,q2,qsurf,
139!       change FF 05/2011
140     .       cloudfrac,totalcloudfrac,hice)
141
142
143
144       ierr = NF_OPEN (fichnom, NF_NOWRITE,nid1)
145       IF (ierr.NE.NF_NOERR) THEN
146         write(6,*)' Pb d''ouverture du fichier'//fichnom
147        CALL ABORT
148       ENDIF
149                                               
150      ierr = NF_INQ_VARID (nid1, "controle", varid)
151      IF (ierr .NE. NF_NOERR) THEN
152       PRINT*, "start2archive: Le champ <controle> est absent"
153       CALL abort
154      ENDIF
155#ifdef NC_DOUBLE
156       ierr = NF_GET_VAR_DOUBLE(nid1, varid, tab_cntrl_fi)
157#else
158      ierr = NF_GET_VAR_REAL(nid1, varid, tab_cntrl_fi)
159#endif
160       IF (ierr .NE. NF_NOERR) THEN
161          PRINT*, "start2archive: Lecture echoue pour <controle>"
162          CALL abort
163       ENDIF
164
165      ierr = NF_CLOSE(nid1)
166
167c-----------------------------------------------------------------------
168c Controle de la synchro
169c-----------------------------------------------------------------------
170!mars a voir      if ((day_ini_fi.ne.day_ini).or.(abs(timefi-timedyn).gt.1.e-10))
171      if ((day_ini_fi.ne.day_ini))
172     &  stop ' Probleme de Synchro entre start et startfi !!!'
173
174
175c *****************************************************************
176c    Option : Reinitialisation des dates dans la premieres annees :
177       do while (day_ini.ge.669)
178          day_ini=day_ini-669
179       enddo
180c *****************************************************************
181
182c-----------------------------------------------------------------------
183c   Initialisations
184c-----------------------------------------------------------------------
185
186      CALL defrun_new(99, .FALSE. )
187      call iniconst
188      call inigeom
189      call inifilr
190      CALL pression(ip1jmp1, ap, bp, ps, p3d)
191      call exner_hyb(ip1jmp1, ps, p3d, beta, pks, pk, pkf)
192
193c=======================================================================
194c Transformation EN VARIABLE NATURELLE / GRILLE SCALAIRE si necessaire
195c=======================================================================
196c  Les variables modeles dependent de la resolution. Il faut donc
197c  eliminer les facteurs responsables de cette dependance
198c  (pour utiliser newstart)
199c=======================================================================
200
201c-----------------------------------------------------------------------
202c Vent   (depend de la resolution horizontale)
203c-----------------------------------------------------------------------
204c
205c ucov --> un  et  vcov --> vn
206c un --> us  et   vn --> vs
207c
208c-----------------------------------------------------------------------
209
210      call covnat(llm,ucov, vcov, un, vn)
211      call wind_scal(un,vn,us,vs)
212
213c-----------------------------------------------------------------------
214c Temperature  (depend de la resolution verticale => de "sigma.def")
215c-----------------------------------------------------------------------
216c
217c h --> T
218c
219c-----------------------------------------------------------------------
220
221      DO l=1,llm
222         DO ij=1,ip1jmp1
223            T(ij,l)=teta(ij,l)*pk(ij,l)/cpp !mars deduit de l'equation dans newstart
224         ENDDO
225      ENDDO
226
227c-----------------------------------------------------------------------
228c Variable physique
229c-----------------------------------------------------------------------
230c
231c tsurf --> tsurfS
232c co2ice --> co2iceS
233c tsoil --> tsoilS
234c emis --> emisS
235c q2 --> q2S
236c qsurf --> qsurfS
237c
238c-----------------------------------------------------------------------
239
240      call gr_fi_dyn(1,ngridmx,iip1,jjp1,tsurf,tsurfS)
241!      call gr_fi_dyn(1,ngridmx,iip1,jjp1,co2ice,co2iceS)
242      call gr_fi_dyn(nsoilmx,ngridmx,iip1,jjp1,tsoil,tsoilS)
243      ! Note: thermal inertia "inertiedat" is in comsoil.h
244      call gr_fi_dyn(nsoilmx,ngridmx,iip1,jjp1,inertiedat,ithS)
245      call gr_fi_dyn(1,ngridmx,iip1,jjp1,emis,emisS)
246      call gr_fi_dyn(llm+1,ngridmx,iip1,jjp1,q2,q2S)
247      call gr_fi_dyn(nqmx,ngridmx,iip1,jjp1,qsurf,qsurfS)
248      call gr_fi_dyn(llm,ngridmx,iip1,jjp1,cloudfrac,cloudfracS)
249      call gr_fi_dyn(1,ngridmx,iip1,jjp1,hice,hiceS)
250      call gr_fi_dyn(1,ngridmx,iip1,jjp1,totalcloudfrac,totalcloudfracS)
251
252c=======================================================================
253c Info pour controler
254c=======================================================================
255
256      ptotal =  0.
257      co2icetotal = 0.
258      DO j=1,jjp1
259         DO i=1,iim
260           ptotal=ptotal+aire(i+(iim+1)*(j-1))*ps(i+(iim+1)*(j-1))/g
261!           co2icetotal = co2icetotal +
262!     &            co2iceS(i+(iim+1)*(j-1))*aire(i+(iim+1)*(j-1))
263         ENDDO
264      ENDDO
265      write(*,*)'Ancienne grille : masse de l''atm :',ptotal
266!      write(*,*)'Ancienne grille : masse de la glace CO2 :',co2icetotal
267
268c-----------------------------------------------------------------------
269c Passage de "ptotal" et "co2icetotal" par tab_cntrl_fi
270c-----------------------------------------------------------------------
271
272      tab_cntrl_fi(49) = ptotal
273      tab_cntrl_fi(50) = co2icetotal
274
275c=======================================================================
276c Ecriture dans le fichier  "start_archive"
277c=======================================================================
278
279c-----------------------------------------------------------------------
280c Ouverture de "start_archive"
281c-----------------------------------------------------------------------
282
283      ierr = NF_OPEN ('start_archive.nc', NF_WRITE,nid)
284 
285c-----------------------------------------------------------------------
286c  si "start_archive" n'existe pas:
287c    1_ ouverture
288c    2_ creation de l'entete dynamique ("ini_archive")
289c-----------------------------------------------------------------------
290c ini_archive:
291c On met dans l'entete le tab_cntrl dynamique (1 a 16)
292c  On y ajoute les valeurs du tab_cntrl_fi (a partir de 51)
293c  En plus les deux valeurs ptotal et co2icetotal (99 et 100)
294c-----------------------------------------------------------------------
295
296      if (ierr.ne.NF_NOERR) then
297         write(*,*)'OK, Could not open file "start_archive.nc"'
298         write(*,*)'So let s create a new "start_archive"'
299         ierr = NF_CREATE('start_archive.nc', NF_CLOBBER, nid)
300         call ini_archive(nid,day_ini,phis,ithS,tab_cntrl_fi)
301      endif
302
303c-----------------------------------------------------------------------
304c Ecriture de la coordonnee temps (date en jours)
305c-----------------------------------------------------------------------
306
307      date = day_ini
308      ierr= NF_INQ_VARID(nid,"Time",varid)
309      ierr= NF_INQ_DIMID(nid,"Time",dimid)
310      ierr= NF_INQ_DIMLEN(nid,dimid,timelen)
311      ntime=timelen+1
312
313      write(*,*) "******************"
314      write(*,*) "ntime",ntime
315      write(*,*) "******************"
316#ifdef NC_DOUBLE
317      ierr= NF_PUT_VARA_DOUBLE(nid,varid,ntime,1,date)
318#else
319      ierr= NF_PUT_VARA_REAL(nid,varid,ntime,1,date)
320#endif
321      if (ierr.ne.NF_NOERR) then
322         write(*,*) "time matter ",NF_STRERROR(ierr)
323         stop
324      endif
325
326c-----------------------------------------------------------------------
327c Ecriture des champs  (co2ice,emis,ps,Tsurf,T,u,v,q2,q,qsurf)
328c-----------------------------------------------------------------------
329c ATTENTION: q2 a une couche de plus!!!!
330c    Pour creer un fichier netcdf lisible par grads,
331c    On passe donc une des couches de q2 a part
332c    comme une variable 2D (la couche au sol: "q2surf")
333c    Les lmm autres couches sont nommees "q2atm" (3D)
334c-----------------------------------------------------------------------
335
336!      call write_archive(nid,ntime,'co2ice','couche de glace co2',
337!     &  'kg/m2',2,co2iceS)
338      call write_archive(nid,ntime,'emis','grd emis',' ',2,emisS)
339      call write_archive(nid,ntime,'ps','Psurf','Pa',2,ps)
340      call write_archive(nid,ntime,'tsurf','surf T','K',2,tsurfS)
341      call write_archive(nid,ntime,'temp','temperature','K',3,t)
342      call write_archive(nid,ntime,'u','Vent zonal','m.s-1',3,us)
343      call write_archive(nid,ntime,'v','Vent merid','m.s-1',3,vs)
344      call write_archive(nid,ntime,'q2surf','wind variance','m2.s-2',2,
345     .              q2S)
346      call write_archive(nid,ntime,'q2atm','wind variance','m2.s-2',3,
347     .              q2S(1,2))
348
349c-----------------------------------------------------------------------
350c Ecriture du champs  q  ( q[1,nqmx] )
351c-----------------------------------------------------------------------
352      do iq=1,nqmx
353        call write_archive(nid,ntime,tnom(iq),'tracer','kg/kg',
354     &         3,q(1,1,iq))
355      end do
356c-----------------------------------------------------------------------
357c Ecriture du champs  qsurf  ( qsurf[1,nqmx] )
358c-----------------------------------------------------------------------
359      do iq=1,nqmx
360        txt=trim(tnom(iq))//"_surf"
361        call write_archive(nid,ntime,txt,'Tracer on surface',
362     &  'kg.m-2',2,qsurfS(1,iq))
363      enddo
364
365
366c-----------------------------------------------------------------------
367c Ecriture du champs  tsoil  ( Tg[1,10] )
368c-----------------------------------------------------------------------
369c "tsoil" Temperature au sol definie dans 10 couches dans le sol
370c   Les 10 couches sont lues comme 10 champs
371c  nommees Tg[1,10]
372
373c      do isoil=1,nsoilmx
374c       write(str2,'(i2.2)') isoil
375c       call write_archive(nid,ntime,'Tg'//str2,'Ground Temperature ',
376c     .   'K',2,tsoilS(1,isoil))
377c      enddo
378
379! Write soil temperatures tsoil
380      call write_archive(nid,ntime,'tsoil','Soil temperature',
381     &     'K',-3,tsoilS)
382
383! Write soil thermal inertia
384      call write_archive(nid,ntime,'inertiedat',
385     &     'Soil thermal inertia',
386     &     'J.s-1/2.m-2.K-1',-3,ithS)
387
388! Write (0D) volumetric heat capacity (stored in comsoil.h)
389!      call write_archive(nid,ntime,'volcapa',
390!     &     'Soil volumetric heat capacity',
391!     &     'J.m-3.K-1',0,volcapa)
392! Note: no need to write volcapa, it is stored in "controle" table
393
394c-----------------------------------------------------------------------
395c Ecriture du champs  cloudfrac,hice,totalcloudfrac
396c-----------------------------------------------------------------------
397      call write_archive(nid,ntime,'hice',
398     &         'Height of oceanic ice','m',2,hiceS)
399      call write_archive(nid,ntime,'totalcloudfrac',
400     &        'Total cloud Fraction','',2,totalcloudfracS)
401      call write_archive(nid,ntime,'cloudfrac'
402     &        ,'Cloud fraction','',3,cloudfracS)
403c-----------------------------------------------------------------------
404c Fin
405c-----------------------------------------------------------------------
406      ierr=NF_CLOSE(nid)
407
408      end
Note: See TracBrowser for help on using the repository browser.