LMDZ
guide_mod.F90
Go to the documentation of this file.
1 !
2 ! $Id$
3 !
4 MODULE guide_mod
5 
6 !=======================================================================
7 ! Auteur: F.Hourdin
8 ! F. Codron 01/09
9 !=======================================================================
10 
11  USE getparam
12  USE write_field
13  use netcdf, only: nf90_nowrite, nf90_open, nf90_inq_varid, nf90_close
14  use pres2lev_mod
15 
16  IMPLICIT NONE
17 
18 ! ---------------------------------------------
19 ! Declarations des cles logiques et parametres
20 ! ---------------------------------------------
21  INTEGER, PRIVATE, SAVE :: iguide_read,iguide_int,iguide_sav
22  INTEGER, PRIVATE, SAVE :: nlevnc
23  LOGICAL, PRIVATE, SAVE :: guide_u,guide_v,guide_t,guide_q,guide_p
24  LOGICAL, PRIVATE, SAVE :: guide_hr,guide_teta
25  LOGICAL, PRIVATE, SAVE :: guide_bl,guide_reg,guide_add,gamma4,guide_zon
26  LOGICAL, PRIVATE, SAVE :: guide_modele,invert_p,invert_y,ini_anal
27  LOGICAL, PRIVATE, SAVE :: guide_2d,guide_sav
28 
29  REAL, PRIVATE, SAVE :: tau_min_u,tau_max_u
30  REAL, PRIVATE, SAVE :: tau_min_v,tau_max_v
31  REAL, PRIVATE, SAVE :: tau_min_t,tau_max_t
32  REAL, PRIVATE, SAVE :: tau_min_q,tau_max_q
33  REAL, PRIVATE, SAVE :: tau_min_p,tau_max_p
34 
35  REAL, PRIVATE, SAVE :: lat_min_g,lat_max_g
36  REAL, PRIVATE, SAVE :: lon_min_g,lon_max_g
37  REAL, PRIVATE, SAVE :: tau_lon,tau_lat
38 
39  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: alpha_u,alpha_v
40  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: alpha_t,alpha_q
41  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: alpha_p,alpha_pcor
42 
43 ! ---------------------------------------------
44 ! Variables de guidage
45 ! ---------------------------------------------
46 ! Variables des fichiers de guidage
47  REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: unat1,unat2
48  REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: vnat1,vnat2
49  REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: tnat1,tnat2
50  REAL, ALLOCATABLE, DIMENSION(:,:,:), PRIVATE, SAVE :: qnat1,qnat2
51  REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: psnat1,psnat2
52  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: apnc,bpnc
53 ! Variables aux dimensions du modele
54  REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: ugui1,ugui2
55  REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: vgui1,vgui2
56  REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: tgui1,tgui2
57  REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: qgui1,qgui2
58  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: psgui1,psgui2
59 
60 CONTAINS
61 !=======================================================================
62 
63  SUBROUTINE guide_init
64 
65  USE control_mod
66 
67  IMPLICIT NONE
68 
69  include "dimensions.h"
70  include "paramet.h"
71  include "netcdf.inc"
72 
73  ! For grossismx:
74  include "serre.h"
75 
76  INTEGER :: error,ncidpl,rid,rcod
77  CHARACTER (len = 80) :: abort_message
78  CHARACTER (len = 20) :: modname = 'guide_init'
79 
80 ! ---------------------------------------------
81 ! Lecture des parametres:
82 ! ---------------------------------------------
83  call ini_getparam("nudging_parameters_out.txt")
84 ! Variables guidees
85  CALL getpar('guide_u',.true.,guide_u,'guidage de u')
86  CALL getpar('guide_v',.true.,guide_v,'guidage de v')
87  CALL getpar('guide_T',.true.,guide_t,'guidage de T')
88  CALL getpar('guide_P',.true.,guide_p,'guidage de P')
89  CALL getpar('guide_Q',.true.,guide_q,'guidage de Q')
90  CALL getpar('guide_hr',.true.,guide_hr,'guidage de Q par H.R')
91  CALL getpar('guide_teta',.false.,guide_teta,'guidage de T par Teta')
92 
93  CALL getpar('guide_add',.false.,guide_add,'for�age constant?')
94  CALL getpar('guide_zon',.false.,guide_zon,'guidage moy zonale')
95  if (guide_zon .and. abs(grossismx - 1.) > 0.01) &
96  call abort_gcm("guide_init", &
97  "zonal nudging requires grid regular in longitude", 1)
98 
99 ! Constantes de rappel. Unite : fraction de jour
100  CALL getpar('tau_min_u',0.02,tau_min_u,'Cste de rappel min, u')
101  CALL getpar('tau_max_u', 10.,tau_max_u,'Cste de rappel max, u')
102  CALL getpar('tau_min_v',0.02,tau_min_v,'Cste de rappel min, v')
103  CALL getpar('tau_max_v', 10.,tau_max_v,'Cste de rappel max, v')
104  CALL getpar('tau_min_T',0.02,tau_min_t,'Cste de rappel min, T')
105  CALL getpar('tau_max_T', 10.,tau_max_t,'Cste de rappel max, T')
106  CALL getpar('tau_min_Q',0.02,tau_min_q,'Cste de rappel min, Q')
107  CALL getpar('tau_max_Q', 10.,tau_max_q,'Cste de rappel max, Q')
108  CALL getpar('tau_min_P',0.02,tau_min_p,'Cste de rappel min, P')
109  CALL getpar('tau_max_P', 10.,tau_max_p,'Cste de rappel max, P')
110  CALL getpar('gamma4',.false.,gamma4,'Zone sans rappel elargie')
111  CALL getpar('guide_BL',.true.,guide_bl,'guidage dans C.Lim')
112 
113 ! Sauvegarde du for�age
114  CALL getpar('guide_sav',.false.,guide_sav,'sauvegarde guidage')
115  CALL getpar('iguide_sav',4,iguide_sav,'freq. sauvegarde guidage')
116  ! frequences f>0: fx/jour; f<0: tous les f jours; f=0: 1 seule fois.
117  IF (iguide_sav.GT.0) THEN
118  iguide_sav=day_step/iguide_sav
119  ELSE if (iguide_sav == 0) then
120  iguide_sav = huge(0)
121  ELSE
122  iguide_sav=day_step*iguide_sav
123  ENDIF
124 
125 ! Guidage regional seulement (sinon constant ou suivant le zoom)
126  CALL getpar('guide_reg',.false.,guide_reg,'guidage regional')
127  CALL getpar('lat_min_g',-90.,lat_min_g,'Latitude mini guidage ')
128  CALL getpar('lat_max_g', 90.,lat_max_g,'Latitude maxi guidage ')
129  CALL getpar('lon_min_g',-180.,lon_min_g,'longitude mini guidage ')
130  CALL getpar('lon_max_g', 180.,lon_max_g,'longitude maxi guidage ')
131  CALL getpar('tau_lat', 5.,tau_lat,'raideur lat guide regional ')
132  CALL getpar('tau_lon', 5.,tau_lon,'raideur lon guide regional ')
133 
134 ! Parametres pour lecture des fichiers
135  CALL getpar('iguide_read',4,iguide_read,'freq. lecture guidage')
136  CALL getpar('iguide_int',4,iguide_int,'freq. interpolation vert')
137  IF (iguide_int.EQ.0) THEN
138  iguide_int=1
139  ELSEIF (iguide_int.GT.0) THEN
140  iguide_int=day_step/iguide_int
141  ELSE
142  iguide_int=day_step*iguide_int
143  ENDIF
144  CALL getpar('guide_modele',.false.,guide_modele,'guidage niveaux modele')
145  CALL getpar('ini_anal',.false.,ini_anal,'Etat initial = analyse')
146  CALL getpar('guide_invertp',.true.,invert_p,'niveaux p inverses')
147  CALL getpar('guide_inverty',.true.,invert_y,'inversion N-S')
148  CALL getpar('guide_2D',.false.,guide_2d,'fichier guidage lat-P')
149 
150  call fin_getparam
151 
152 ! ---------------------------------------------
153 ! Determination du nombre de niveaux verticaux
154 ! des fichiers guidage
155 ! ---------------------------------------------
156  ncidpl=-99
157  if (guide_modele) then
158  if (ncidpl.eq.-99) then
159  rcod=nf90_open('apbp.nc',nf90_nowrite, ncidpl)
160  if (rcod.NE.nf_noerr) THEN
161  CALL abort_gcm(modname, &
162  'Guide: probleme -> pas de fichier apbp.nc',1)
163  endif
164  endif
165  else
166  if (guide_u) then
167  if (ncidpl.eq.-99) then
168  rcod=nf90_open('u.nc',nf90_nowrite,ncidpl)
169  if (rcod.NE.nf_noerr) THEN
170  CALL abort_gcm(modname, &
171  'Guide: probleme -> pas de fichier u.nc',1)
172  endif
173  endif
174  elseif (guide_v) then
175  if (ncidpl.eq.-99) then
176  rcod=nf90_open('v.nc',nf90_nowrite,ncidpl)
177  if (rcod.NE.nf_noerr) THEN
178  CALL abort_gcm(modname, &
179  'Guide: probleme -> pas de fichier v.nc',1)
180  endif
181  endif
182  elseif (guide_t) then
183  if (ncidpl.eq.-99) then
184  rcod=nf90_open('T.nc',nf90_nowrite,ncidpl)
185  if (rcod.NE.nf_noerr) THEN
186  CALL abort_gcm(modname, &
187  'Guide: probleme -> pas de fichier T.nc',1)
188  endif
189  endif
190  elseif (guide_q) then
191  if (ncidpl.eq.-99) then
192  rcod=nf90_open('hur.nc',nf90_nowrite, ncidpl)
193  if (rcod.NE.nf_noerr) THEN
194  CALL abort_gcm(modname, &
195  'Guide: probleme -> pas de fichier hur.nc',1)
196  endif
197  endif
198  endif
199  endif
200  error=nf_inq_dimid(ncidpl,'LEVEL',rid)
201  IF (error.NE.nf_noerr) error=nf_inq_dimid(ncidpl,'PRESSURE',rid)
202  IF (error.NE.nf_noerr) THEN
203  CALL abort_gcm(modname,'Guide: probleme lecture niveaux pression',1)
204  ENDIF
205  error=nf_inq_dimlen(ncidpl,rid,nlevnc)
206  print *,'Guide: nombre niveaux vert. nlevnc', nlevnc
207  rcod = nf90_close(ncidpl)
208 
209 ! ---------------------------------------------
210 ! Allocation des variables
211 ! ---------------------------------------------
212  abort_message='pb in allocation guide'
213 
214  ALLOCATE(apnc(nlevnc), stat = error)
215  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
216  ALLOCATE(bpnc(nlevnc), stat = error)
217  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
218  apnc=0.;bpnc=0.
219 
220  ALLOCATE(alpha_pcor(llm), stat = error)
221  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
222  ALLOCATE(alpha_u(ip1jmp1), stat = error)
223  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
224  ALLOCATE(alpha_v(ip1jm), stat = error)
225  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
226  ALLOCATE(alpha_t(ip1jmp1), stat = error)
227  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
228  ALLOCATE(alpha_q(ip1jmp1), stat = error)
229  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
230  ALLOCATE(alpha_p(ip1jmp1), stat = error)
231  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
232  alpha_u=0.;alpha_v=0;alpha_t=0;alpha_q=0;alpha_p=0
233 
234  IF (guide_u) THEN
235  ALLOCATE(unat1(iip1,jjp1,nlevnc), stat = error)
236  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
237  ALLOCATE(ugui1(ip1jmp1,llm), stat = error)
238  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
239  ALLOCATE(unat2(iip1,jjp1,nlevnc), stat = error)
240  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
241  ALLOCATE(ugui2(ip1jmp1,llm), stat = error)
242  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
243  unat1=0.;unat2=0.;ugui1=0.;ugui2=0.
244  ENDIF
245 
246  IF (guide_t) THEN
247  ALLOCATE(tnat1(iip1,jjp1,nlevnc), stat = error)
248  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
249  ALLOCATE(tgui1(ip1jmp1,llm), stat = error)
250  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
251  ALLOCATE(tnat2(iip1,jjp1,nlevnc), stat = error)
252  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
253  ALLOCATE(tgui2(ip1jmp1,llm), stat = error)
254  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
255  tnat1=0.;tnat2=0.;tgui1=0.;tgui2=0.
256  ENDIF
257 
258  IF (guide_q) THEN
259  ALLOCATE(qnat1(iip1,jjp1,nlevnc), stat = error)
260  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
261  ALLOCATE(qgui1(ip1jmp1,llm), stat = error)
262  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
263  ALLOCATE(qnat2(iip1,jjp1,nlevnc), stat = error)
264  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
265  ALLOCATE(qgui2(ip1jmp1,llm), stat = error)
266  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
267  qnat1=0.;qnat2=0.;qgui1=0.;qgui2=0.
268  ENDIF
269 
270  IF (guide_v) THEN
271  ALLOCATE(vnat1(iip1,jjm,nlevnc), stat = error)
272  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
273  ALLOCATE(vgui1(ip1jm,llm), stat = error)
274  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
275  ALLOCATE(vnat2(iip1,jjm,nlevnc), stat = error)
276  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
277  ALLOCATE(vgui2(ip1jm,llm), stat = error)
278  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
279  vnat1=0.;vnat2=0.;vgui1=0.;vgui2=0.
280  ENDIF
281 
282  IF (guide_p.OR.guide_modele) THEN
283  ALLOCATE(psnat1(iip1,jjp1), stat = error)
284  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
285  ALLOCATE(psnat2(iip1,jjp1), stat = error)
286  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
287  psnat1=0.;psnat2=0.;
288  ENDIF
289  IF (guide_p) THEN
290  ALLOCATE(psgui2(ip1jmp1), stat = error)
291  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
292  ALLOCATE(psgui1(ip1jmp1), stat = error)
293  IF (error /= 0) CALL abort_gcm(modname,abort_message,1)
294  psgui1=0.;psgui2=0.
295  ENDIF
296 
297 ! ---------------------------------------------
298 ! Lecture du premier etat de guidage.
299 ! ---------------------------------------------
300  IF (guide_2d) THEN
301  CALL guide_read2d(1)
302  ELSE
303  CALL guide_read(1)
304  ENDIF
305  IF (guide_v) vnat1=vnat2
306  IF (guide_u) unat1=unat2
307  IF (guide_t) tnat1=tnat2
308  IF (guide_q) qnat1=qnat2
309  IF (guide_p.OR.guide_modele) psnat1=psnat2
310 
311  END SUBROUTINE guide_init
312 
313 !=======================================================================
314  SUBROUTINE guide_main(itau,ucov,vcov,teta,q,masse,ps)
315 
316  USE control_mod
317 
318  IMPLICIT NONE
319 
320  include "dimensions.h"
321  include "paramet.h"
322  include "comconst.h"
323  include "comvert.h"
324 
325  ! Variables entree
326  INTEGER, INTENT(IN) :: itau !pas de temps
327  REAL, DIMENSION (ip1jmp1,llm), INTENT(INOUT) :: ucov,teta,q,masse
328  REAL, DIMENSION (ip1jm,llm), INTENT(INOUT) :: vcov
329  REAL, DIMENSION (ip1jmp1), INTENT(INOUT) :: ps
330 
331  ! Variables locales
332  LOGICAL, SAVE :: first=.true.
333  LOGICAL :: f_out ! sortie guidage
334  REAL, DIMENSION (ip1jmp1,llm) :: f_add ! var aux: champ de guidage
335  REAL, DIMENSION (ip1jmp1,llm) :: p ! besoin si guide_P
336  ! Compteurs temps:
337  INTEGER, SAVE :: step_rea,count_no_rea,itau_test ! lecture guidage
338  REAL :: ditau, dday_step
339  REAL :: tau,reste ! position entre 2 etats de guidage
340  REAL, SAVE :: factt ! pas de temps en fraction de jour
341 
342  INTEGER :: l
343 
344 !-----------------------------------------------------------------------
345 ! Initialisations au premier passage
346 !-----------------------------------------------------------------------
347  IF (first) THEN
348  first=.false.
349  CALL guide_init
350  itau_test=1001
351  step_rea=1
352  count_no_rea=0
353 ! Calcul des constantes de rappel
354  factt=dtvr*iperiod/daysec
355  call tau2alpha(3,iip1,jjm ,factt,tau_min_v,tau_max_v,alpha_v)
356  call tau2alpha(2,iip1,jjp1,factt,tau_min_u,tau_max_u,alpha_u)
357  call tau2alpha(1,iip1,jjp1,factt,tau_min_t,tau_max_t,alpha_t)
358  call tau2alpha(1,iip1,jjp1,factt,tau_min_p,tau_max_p,alpha_p)
359  call tau2alpha(1,iip1,jjp1,factt,tau_min_q,tau_max_q,alpha_q)
360 ! correction de rappel dans couche limite
361  if (guide_bl) then
362  alpha_pcor(:)=1.
363  else
364  do l=1,llm
365  alpha_pcor(l)=(1.+tanh((0.85-presnivs(l)/preff)/0.05))/2.
366  enddo
367  endif
368 ! ini_anal: etat initial egal au guidage
369  IF (ini_anal) THEN
370  CALL guide_interp(ps,teta)
371  IF (guide_u) ucov=ugui2
372  IF (guide_v) vcov=ugui2
373  IF (guide_t) teta=tgui2
374  IF (guide_q) q=qgui2
375  IF (guide_p) THEN
376  ps=psgui2
377  CALL pression(ip1jmp1,ap,bp,ps,p)
378  CALL massdair(p,masse)
379  ENDIF
380  RETURN
381  ENDIF
382 ! Verification structure guidage
383  IF (guide_u) THEN
384  CALL writefield('unat',unat1)
385  CALL writefield('ucov',reshape(ucov,(/iip1,jjp1,llm/)))
386  ENDIF
387  IF (guide_t) THEN
388  CALL writefield('tnat',tnat1)
389  CALL writefield('teta',reshape(teta,(/iip1,jjp1,llm/)))
390  ENDIF
391 
392  ENDIF !first
393 
394 !-----------------------------------------------------------------------
395 ! Lecture des fichiers de guidage ?
396 !-----------------------------------------------------------------------
397  IF (iguide_read.NE.0) THEN
398  ditau=real(itau)
399  dday_step=real(day_step)
400  IF (iguide_read.LT.0) THEN
401  tau=ditau/dday_step/REAL(iguide_read)
402  ELSE
403  tau=REAL(iguide_read)*ditau/dday_step
404  ENDIF
405  reste=tau-aint(tau)
406  IF (reste.EQ.0.) THEN
407  IF (itau_test.EQ.itau) THEN
408  write(*,*)'deuxieme passage de advreel a itau=',itau
409  stop
410  ELSE
411  IF (guide_v) vnat1=vnat2
412  IF (guide_u) unat1=unat2
413  IF (guide_t) tnat1=tnat2
414  IF (guide_q) qnat1=qnat2
415  IF (guide_p.OR.guide_modele) psnat1=psnat2
416  step_rea=step_rea+1
417  itau_test=itau
418  print*,'Lecture fichiers guidage, pas ',step_rea, &
419  'apres ',count_no_rea,' non lectures'
420  IF (guide_2d) THEN
421  CALL guide_read2d(step_rea)
422  ELSE
423  CALL guide_read(step_rea)
424  ENDIF
425  count_no_rea=0
426  ENDIF
427  ELSE
428  count_no_rea=count_no_rea+1
429 
430  ENDIF
431  ENDIF !iguide_read=0
432 
433 !-----------------------------------------------------------------------
434 ! Interpolation et conversion des champs de guidage
435 !-----------------------------------------------------------------------
436  IF (mod(itau,iguide_int).EQ.0) THEN
437  CALL guide_interp(ps,teta)
438  ENDIF
439 ! Repartition entre 2 etats de guidage
440  IF (iguide_read.NE.0) THEN
441  tau=reste
442  ELSE
443  tau=1.
444  ENDIF
445 
446 !-----------------------------------------------------------------------
447 ! Ajout des champs de guidage
448 !-----------------------------------------------------------------------
449 ! Sauvegarde du guidage?
450  f_out=((mod(itau,iguide_sav).EQ.0).AND.guide_sav)
451  IF (f_out) CALL guide_out("SP",jjp1,1,ps)
452 
453  if (guide_u) then
454  if (guide_add) then
455  f_add=(1.-tau)*ugui1+tau*ugui2
456  else
457  f_add=(1.-tau)*ugui1+tau*ugui2-ucov
458  endif
459  if (guide_zon) CALL guide_zonave(1,jjp1,llm,f_add)
460  CALL guide_addfield(ip1jmp1,llm,f_add,alpha_u)
461  IF (f_out) CALL guide_out("ua",jjp1,llm,(1.-tau)*ugui1+tau*ugui2)
462  IF (f_out) CALL guide_out("u",jjp1,llm,ucov)
463  IF (f_out) CALL guide_out("ucov",jjp1,llm,f_add/factt)
464  ucov=ucov+f_add
465  endif
466 
467  if (guide_t) then
468  if (guide_add) then
469  f_add=(1.-tau)*tgui1+tau*tgui2
470  else
471  f_add=(1.-tau)*tgui1+tau*tgui2-teta
472  endif
473  if (guide_zon) CALL guide_zonave(2,jjp1,llm,f_add)
474  CALL guide_addfield(ip1jmp1,llm,f_add,alpha_t)
475  IF (f_out) CALL guide_out("teta",jjp1,llm,f_add/factt)
476  teta=teta+f_add
477  endif
478 
479  if (guide_p) then
480  if (guide_add) then
481  f_add(1:ip1jmp1,1)=(1.-tau)*psgui1+tau*psgui2
482  else
483  f_add(1:ip1jmp1,1)=(1.-tau)*psgui1+tau*psgui2-ps
484  endif
485  if (guide_zon) CALL guide_zonave(2,jjp1,1,f_add(1:ip1jmp1,1))
486  CALL guide_addfield(ip1jmp1,1,f_add(1:ip1jmp1,1),alpha_p)
487  IF (f_out) CALL guide_out("ps",jjp1,1,f_add(1:ip1jmp1,1)/factt)
488  ps=ps+f_add(1:ip1jmp1,1)
489  CALL pression(ip1jmp1,ap,bp,ps,p)
490  CALL massdair(p,masse)
491  endif
492 
493  if (guide_q) then
494  if (guide_add) then
495  f_add=(1.-tau)*qgui1+tau*qgui2
496  else
497  f_add=(1.-tau)*qgui1+tau*qgui2-q
498  endif
499  if (guide_zon) CALL guide_zonave(2,jjp1,llm,f_add)
500  CALL guide_addfield(ip1jmp1,llm,f_add,alpha_q)
501  IF (f_out) CALL guide_out("q",jjp1,llm,f_add/factt)
502  q=q+f_add
503  endif
504 
505  if (guide_v) then
506  if (guide_add) then
507  f_add(1:ip1jm,:)=(1.-tau)*vgui1+tau*vgui2
508  else
509  f_add(1:ip1jm,:)=(1.-tau)*vgui1+tau*vgui2-vcov
510  endif
511  if (guide_zon) CALL guide_zonave(2,jjm,llm,f_add(1:ip1jm,:))
512  CALL guide_addfield(ip1jm,llm,f_add(1:ip1jm,:),alpha_v)
513  IF (f_out) CALL guide_out("v",jjm,llm,vcov)
514  IF (f_out) CALL guide_out("va",jjm,llm,(1.-tau)*vgui1+tau*vgui2)
515  IF (f_out) CALL guide_out("vcov",jjm,llm,f_add(1:ip1jm,:)/factt)
516  vcov=vcov+f_add(1:ip1jm,:)
517  endif
518  END SUBROUTINE guide_main
519 
520 !=======================================================================
521  SUBROUTINE guide_addfield(hsize,vsize,field,alpha)
522 ! field1=a*field1+alpha*field2
523 
524  IMPLICIT NONE
525 
526  ! input variables
527  INTEGER, INTENT(IN) :: hsize
528  INTEGER, INTENT(IN) :: vsize
529  REAL, DIMENSION(hsize), INTENT(IN) :: alpha
530  REAL, DIMENSION(hsize,vsize), INTENT(INOUT) :: field
531 
532  ! Local variables
533  INTEGER :: l
534 
535  do l=1,vsize
536  field(:,l)=alpha*field(:,l)*alpha_pcor(l)
537  enddo
538 
539  END SUBROUTINE guide_addfield
540 
541 !=======================================================================
542  SUBROUTINE guide_zonave(typ,hsize,vsize,field)
543 
544  IMPLICIT NONE
545 
546  include "dimensions.h"
547  include "paramet.h"
548  include "comgeom.h"
549  include "comconst.h"
550 
551  ! input/output variables
552  INTEGER, INTENT(IN) :: typ
553  INTEGER, INTENT(IN) :: vsize
554  INTEGER, INTENT(IN) :: hsize
555  REAL, DIMENSION(hsize*iip1,vsize), INTENT(INOUT) :: field
556 
557  ! Local variables
558  LOGICAL, SAVE :: first=.true.
559  INTEGER, DIMENSION (2), SAVE :: imin, imax ! averaging domain
560  INTEGER :: i,j,l,ij
561  REAL, DIMENSION (iip1) :: lond ! longitude in Deg.
562  REAL, DIMENSION (hsize,vsize):: fieldm ! zon-averaged field
563 
564  IF (first) THEN
565  first=.false.
566 !Compute domain for averaging
567  lond=rlonu*180./pi
568  imin(1)=1;imax(1)=iip1;
569  imin(2)=1;imax(2)=iip1;
570  IF (guide_reg) THEN
571  DO i=1,iim
572  IF (lond(i).LT.lon_min_g) imin(1)=i
573  IF (lond(i).LE.lon_max_g) imax(1)=i
574  ENDDO
575  lond=rlonv*180./pi
576  DO i=1,iim
577  IF (lond(i).LT.lon_min_g) imin(2)=i
578  IF (lond(i).LE.lon_max_g) imax(2)=i
579  ENDDO
580  ENDIF
581  ENDIF
582 
583  fieldm=0.
584  DO l=1,vsize
585  ! Compute zonal average
586  DO j=1,hsize
587  DO i=imin(typ),imax(typ)
588  ij=(j-1)*iip1+i
589  fieldm(j,l)=fieldm(j,l)+field(ij,l)
590  ENDDO
591  ENDDO
592  fieldm(:,l)=fieldm(:,l)/REAL(imax(typ)-imin(typ)+1)
593  ! Compute forcing
594  DO j=1,hsize
595  DO i=1,iip1
596  ij=(j-1)*iip1+i
597  field(ij,l)=fieldm(j,l)
598  ENDDO
599  ENDDO
600  ENDDO
601 
602  END SUBROUTINE guide_zonave
603 
604 !=======================================================================
605  SUBROUTINE guide_interp(psi,teta)
606 
607  use exner_hyb_m, only: exner_hyb
608  use exner_milieu_m, only: exner_milieu
609  IMPLICIT NONE
610 
611  include "dimensions.h"
612  include "paramet.h"
613  include "comvert.h"
614  include "comgeom2.h"
615  include "comconst.h"
616 
617  REAL, DIMENSION (iip1,jjp1), INTENT(IN) :: psi ! Psol gcm
618  REAL, DIMENSION (iip1,jjp1,llm), INTENT(IN) :: teta ! Temp. Pot. gcm
619 
620  LOGICAL, SAVE :: first=.true.
621  ! Variables pour niveaux pression:
622  REAL, DIMENSION (iip1,jjp1,nlevnc) :: plnc1,plnc2 !niveaux pression guidage
623  REAL, DIMENSION (iip1,jjp1,llm) :: plunc,plsnc !niveaux pression modele
624  REAL, DIMENSION (iip1,jjm,llm) :: plvnc !niveaux pression modele
625  REAL, DIMENSION (iip1,jjp1,llmp1) :: p ! pression intercouches
626  REAL, DIMENSION (iip1,jjp1,llm) :: pls, pext ! var intermediaire
627  REAL, DIMENSION (iip1,jjp1,llm) :: pbarx
628  REAL, DIMENSION (iip1,jjm,llm) :: pbary
629  ! Variables pour fonction Exner (P milieu couche)
630  REAL, DIMENSION (iip1,jjp1,llm) :: pk
631  REAL, DIMENSION (iip1,jjp1) :: pks
632  REAL :: prefkap,unskap
633  ! Pression de vapeur saturante
634  REAL, DIMENSION (ip1jmp1,llm) :: qsat
635  !Variables intermediaires interpolation
636  REAL, DIMENSION (iip1,jjp1,llm) :: zu1,zu2
637  REAL, DIMENSION (iip1,jjm,llm) :: zv1,zv2
638 
639  INTEGER :: i,j,l,ij
640 
641  print *,'Guide: conversion variables guidage'
642 ! -----------------------------------------------------------------
643 ! Calcul des niveaux de pression champs guidage
644 ! -----------------------------------------------------------------
645 if (guide_modele) then
646  do i=1,iip1
647  do j=1,jjp1
648  do l=1,nlevnc
649  plnc2(i,j,l)=apnc(l)+bpnc(l)*psnat2(i,j)
650  plnc1(i,j,l)=apnc(l)+bpnc(l)*psnat1(i,j)
651  enddo
652  enddo
653  enddo
654 else
655  do i=1,iip1
656  do j=1,jjp1
657  do l=1,nlevnc
658  plnc2(i,j,l)=apnc(l)
659  plnc1(i,j,l)=apnc(l)
660  enddo
661  enddo
662  enddo
663 
664 endif
665  if (first) then
666  first=.false.
667  print*,'Guide: verification ordre niveaux verticaux'
668  print*,'LMDZ :'
669  do l=1,llm
670  print*,'PL(',l,')=',(ap(l)+ap(l+1))/2. &
671  +psi(1,jjp1)*(bp(l)+bp(l+1))/2.
672  enddo
673  print*,'Fichiers guidage'
674  do l=1,nlevnc
675  print*,'PL(',l,')=',plnc2(1,1,l)
676  enddo
677  print *,'inversion de l''ordre: invert_p=',invert_p
678  if (guide_u) then
679  do l=1,nlevnc
680  print*,'U(',l,')=',unat2(1,1,l)
681  enddo
682  endif
683  if (guide_t) then
684  do l=1,nlevnc
685  print*,'T(',l,')=',tnat2(1,1,l)
686  enddo
687  endif
688  endif
689 
690 ! -----------------------------------------------------------------
691 ! Calcul niveaux pression modele
692 ! -----------------------------------------------------------------
693  CALL pression( ip1jmp1, ap, bp, psi, p )
694  if (pressure_exner) then
695  CALL exner_hyb(ip1jmp1,psi,p,pks,pk)
696  else
697  CALL exner_milieu(ip1jmp1,psi,p,pks,pk)
698  endif
699 ! .... Calcul de pls , pression au milieu des couches ,en Pascals
700  unskap=1./kappa
701  prefkap = preff ** kappa
702  DO l = 1, llm
703  DO j=1,jjp1
704  DO i =1, iip1
705  pls(i,j,l) = preff * ( pk(i,j,l)/cpp) ** unskap
706  ENDDO
707  ENDDO
708  ENDDO
709 
710 ! calcul des pressions pour les grilles u et v
711  do l=1,llm
712  do j=1,jjp1
713  do i=1,iip1
714  pext(i,j,l)=pls(i,j,l)*aire(i,j)
715  enddo
716  enddo
717  enddo
718  call massbar(pext, pbarx, pbary )
719  do l=1,llm
720  do j=1,jjp1
721  do i=1,iip1
722  plunc(i,j,l)=pbarx(i,j,l)/aireu(i,j)
723  plsnc(i,j,l)=pls(i,j,l)
724  enddo
725  enddo
726  enddo
727  do l=1,llm
728  do j=1,jjm
729  do i=1,iip1
730  plvnc(i,j,l)=pbary(i,j,l)/airev(i,j)
731  enddo
732  enddo
733  enddo
734 
735 ! -----------------------------------------------------------------
736 ! Interpolation champs guidage sur niveaux modele (+inversion N/S)
737 ! Conversion en variables gcm (ucov, vcov...)
738 ! -----------------------------------------------------------------
739  if (guide_p) then
740  do j=1,jjp1
741  do i=1,iim
742  ij=(j-1)*iip1+i
743  psgui1(ij)=psnat1(i,j)
744  psgui2(ij)=psnat2(i,j)
745  enddo
746  psgui1(iip1*j)=psnat1(1,j)
747  psgui2(iip1*j)=psnat2(1,j)
748  enddo
749  endif
750 
751  IF (guide_u) THEN
752  CALL pres2lev(unat1,zu1,nlevnc,llm,plnc1,plunc,iip1,jjp1,invert_p)
753  CALL pres2lev(unat2,zu2,nlevnc,llm,plnc2,plunc,iip1,jjp1,invert_p)
754  do l=1,llm
755  do j=1,jjp1
756  do i=1,iim
757  ij=(j-1)*iip1+i
758  ugui1(ij,l)=zu1(i,j,l)*cu(i,j)
759  ugui2(ij,l)=zu2(i,j,l)*cu(i,j)
760  enddo
761  ugui1(j*iip1,l)=ugui1((j-1)*iip1+1,l)
762  ugui2(j*iip1,l)=ugui2((j-1)*iip1+1,l)
763  enddo
764  do i=1,iip1
765  ugui1(i,l)=0.
766  ugui1(ip1jm+i,l)=0.
767  ugui2(i,l)=0.
768  ugui2(ip1jm+i,l)=0.
769  enddo
770  enddo
771  ENDIF
772 
773  IF (guide_t) THEN
774  CALL pres2lev(tnat1,zu1,nlevnc,llm,plnc1,plsnc,iip1,jjp1,invert_p)
775  CALL pres2lev(tnat2,zu2,nlevnc,llm,plnc2,plsnc,iip1,jjp1,invert_p)
776  do l=1,llm
777  do j=1,jjp1
778  IF (guide_teta) THEN
779  do i=1,iim
780  ij=(j-1)*iip1+i
781  tgui1(ij,l)=zu1(i,j,l)
782  tgui2(ij,l)=zu2(i,j,l)
783  enddo
784  ELSE
785  do i=1,iim
786  ij=(j-1)*iip1+i
787  tgui1(ij,l)=zu1(i,j,l)*cpp/pk(i,j,l)
788  tgui2(ij,l)=zu2(i,j,l)*cpp/pk(i,j,l)
789  enddo
790  ENDIF
791  tgui1(j*iip1,l)=tgui1((j-1)*iip1+1,l)
792  tgui2(j*iip1,l)=tgui2((j-1)*iip1+1,l)
793  enddo
794  do i=1,iip1
795  tgui1(i,l)=tgui1(1,l)
796  tgui1(ip1jm+i,l)=tgui1(ip1jm+1,l)
797  tgui2(i,l)=tgui2(1,l)
798  tgui2(ip1jm+i,l)=tgui2(ip1jm+1,l)
799  enddo
800  enddo
801  ENDIF
802 
803  IF (guide_v) THEN
804 
805  CALL pres2lev(vnat1,zv1,nlevnc,llm,plnc1(:,:jjm,:),plvnc,iip1,jjm,invert_p)
806  CALL pres2lev(vnat2,zv2,nlevnc,llm,plnc2(:,:jjm,:),plvnc,iip1,jjm,invert_p)
807 
808  do l=1,llm
809  do j=1,jjm
810  do i=1,iim
811  ij=(j-1)*iip1+i
812  vgui1(ij,l)=zv1(i,j,l)*cv(i,j)
813  vgui2(ij,l)=zv2(i,j,l)*cv(i,j)
814  enddo
815  vgui1(j*iip1,l)=vgui1((j-1)*iip1+1,l)
816  vgui2(j*iip1,l)=vgui2((j-1)*iip1+1,l)
817  enddo
818  enddo
819  ENDIF
820 
821  IF (guide_q) THEN
822  ! On suppose qu'on a la bonne variable dans le fichier de guidage:
823  ! Hum.Rel si guide_hr, Hum.Spec. sinon.
824  CALL pres2lev(qnat1,zu1,nlevnc,llm,plnc1,plsnc,iip1,jjp1,invert_p)
825  CALL pres2lev(qnat2,zu2,nlevnc,llm,plnc2,plsnc,iip1,jjp1,invert_p)
826  do l=1,llm
827  do j=1,jjp1
828  do i=1,iim
829  ij=(j-1)*iip1+i
830  qgui1(ij,l)=zu1(i,j,l)
831  qgui2(ij,l)=zu2(i,j,l)
832  enddo
833  qgui1(j*iip1,l)=qgui1((j-1)*iip1+1,l)
834  qgui2(j*iip1,l)=qgui2((j-1)*iip1+1,l)
835  enddo
836  do i=1,iip1
837  qgui1(i,l)=qgui1(1,l)
838  qgui1(ip1jm+i,l)=qgui1(ip1jm+1,l)
839  qgui2(i,l)=qgui2(1,l)
840  qgui2(ip1jm+i,l)=qgui2(ip1jm+1,l)
841  enddo
842  enddo
843  IF (guide_hr) THEN
844  CALL q_sat(iip1*jjp1*llm,teta*pk/cpp,plsnc,qsat)
845  qgui1=qgui1*qsat*0.01 !hum. rel. en %
846  qgui2=qgui2*qsat*0.01
847  ENDIF
848  ENDIF
849 
850  END SUBROUTINE guide_interp
851 
852 !=======================================================================
853  SUBROUTINE tau2alpha(typ,pim,pjm,factt,taumin,taumax,alpha)
854 
855 ! Calcul des constantes de rappel alpha (=1/tau)
856 
857  implicit none
858 
859  include "dimensions.h"
860  include "paramet.h"
861  include "comconst.h"
862  include "comgeom2.h"
863  include "serre.h"
864 
865 ! input arguments :
866  INTEGER, INTENT(IN) :: typ ! u(2),v(3), ou scalaire(1)
867  INTEGER, INTENT(IN) :: pim,pjm ! dimensions en lat, lon
868  REAL, INTENT(IN) :: factt ! pas de temps en fraction de jour
869  REAL, INTENT(IN) :: taumin,taumax
870 ! output arguments:
871  REAL, DIMENSION(pim,pjm), INTENT(OUT) :: alpha
872 
873 ! local variables:
874  LOGICAL, SAVE :: first=.true.
875  REAL, SAVE :: gamma,dxdy_min,dxdy_max
876  REAL, DIMENSION (iip1,jjp1) :: zdx,zdy
877  REAL, DIMENSION (iip1,jjp1) :: dxdys,dxdyu
878  REAL, DIMENSION (iip1,jjm) :: dxdyv
879  real dxdy_
880  real zlat,zlon
881  real alphamin,alphamax,xi
882  integer i,j,ilon,ilat
883 
884 
885  alphamin=factt/taumax
886  alphamax=factt/taumin
887  IF (guide_reg.OR.guide_add) THEN
888  alpha=alphamax
889 !-----------------------------------------------------------------------
890 ! guide_reg: alpha=alpha_min dans region, 0. sinon.
891 !-----------------------------------------------------------------------
892  IF (guide_reg) THEN
893  do j=1,pjm
894  do i=1,pim
895  if (typ.eq.2) then
896  zlat=rlatu(j)*180./pi
897  zlon=rlonu(i)*180./pi
898  elseif (typ.eq.1) then
899  zlat=rlatu(j)*180./pi
900  zlon=rlonv(i)*180./pi
901  elseif (typ.eq.3) then
902  zlat=rlatv(j)*180./pi
903  zlon=rlonv(i)*180./pi
904  endif
905  alpha(i,j)=alphamax/16.* &
906  (1.+tanh((zlat-lat_min_g)/tau_lat))* &
907  (1.+tanh((lat_max_g-zlat)/tau_lat))* &
908  (1.+tanh((zlon-lon_min_g)/tau_lon))* &
909  (1.+tanh((lon_max_g-zlon)/tau_lon))
910  enddo
911  enddo
912  ENDIF
913  ELSE
914 !-----------------------------------------------------------------------
915 ! Sinon, alpha varie entre alpha_min et alpha_max suivant le zoom.
916 !-----------------------------------------------------------------------
917 !Calcul de l'aire des mailles
918  do j=2,jjm
919  do i=2,iip1
920  zdx(i,j)=0.5*(cu(i-1,j)+cu(i,j))/cos(rlatu(j))
921  enddo
922  zdx(1,j)=zdx(iip1,j)
923  enddo
924  do j=2,jjm
925  do i=1,iip1
926  zdy(i,j)=0.5*(cv(i,j-1)+cv(i,j))
927  enddo
928  enddo
929  do i=1,iip1
930  zdx(i,1)=zdx(i,2)
931  zdx(i,jjp1)=zdx(i,jjm)
932  zdy(i,1)=zdy(i,2)
933  zdy(i,jjp1)=zdy(i,jjm)
934  enddo
935  do j=1,jjp1
936  do i=1,iip1
937  dxdys(i,j)=sqrt(zdx(i,j)*zdx(i,j)+zdy(i,j)*zdy(i,j))
938  enddo
939  enddo
940  IF (typ.EQ.2) THEN
941  do j=1,jjp1
942  do i=1,iim
943  dxdyu(i,j)=0.5*(dxdys(i,j)+dxdys(i+1,j))
944  enddo
945  dxdyu(iip1,j)=dxdyu(1,j)
946  enddo
947  ENDIF
948  IF (typ.EQ.3) THEN
949  do j=1,jjm
950  do i=1,iip1
951  dxdyv(i,j)=0.5*(dxdys(i,j)+dxdys(i,j+1))
952  enddo
953  enddo
954  ENDIF
955 ! Premier appel: calcul des aires min et max et de gamma.
956  IF (first) THEN
957  first=.false.
958  ! coordonnees du centre du zoom
959  CALL coordij(clon,clat,ilon,ilat)
960  ! aire de la maille au centre du zoom
961  dxdy_min=dxdys(ilon,ilat)
962  ! dxdy maximale de la maille
963  dxdy_max=0.
964  do j=1,jjp1
965  do i=1,iip1
966  dxdy_max=max(dxdy_max,dxdys(i,j))
967  enddo
968  enddo
969  ! Calcul de gamma
970  if (abs(grossismx-1.).lt.0.1.or.abs(grossismy-1.).lt.0.1) then
971  print*,'ATTENTION modele peu zoome'
972  print*,'ATTENTION on prend une constante de guidage cste'
973  gamma=0.
974  else
975  gamma=(dxdy_max-2.*dxdy_min)/(dxdy_max-dxdy_min)
976  print*,'gamma=',gamma
977  if (gamma.lt.1.e-5) then
978  print*,'gamma =',gamma,'<1e-5'
979  stop
980  endif
981  gamma=log(0.5)/log(gamma)
982  if (gamma4) then
983  gamma=min(gamma,4.)
984  endif
985  print*,'gamma=',gamma
986  endif
987  ENDIF !first
988 
989  do j=1,pjm
990  do i=1,pim
991  if (typ.eq.1) then
992  dxdy_=dxdys(i,j)
993  zlat=rlatu(j)*180./pi
994  elseif (typ.eq.2) then
995  dxdy_=dxdyu(i,j)
996  zlat=rlatu(j)*180./pi
997  elseif (typ.eq.3) then
998  dxdy_=dxdyv(i,j)
999  zlat=rlatv(j)*180./pi
1000  endif
1001  if (abs(grossismx-1.).lt.0.1.or.abs(grossismy-1.).lt.0.1) then
1002  ! pour une grille reguliere, xi=xxx**0=1 -> alpha=alphamin
1003  alpha(i,j)=alphamin
1004  else
1005  xi=((dxdy_max-dxdy_)/(dxdy_max-dxdy_min))**gamma
1006  xi=min(xi,1.)
1007  if(lat_min_g.le.zlat .and. zlat.le.lat_max_g) then
1008  alpha(i,j)=xi*alphamin+(1.-xi)*alphamax
1009  else
1010  alpha(i,j)=0.
1011  endif
1012  endif
1013  enddo
1014  enddo
1015  ENDIF ! guide_reg
1016 
1017  if (.not. guide_add) alpha = 1. - exp(- alpha)
1018 
1019  END SUBROUTINE tau2alpha
1020 
1021 !=======================================================================
1022  SUBROUTINE guide_read(timestep)
1023 
1024  IMPLICIT NONE
1025 
1026 #include "netcdf.inc"
1027 #include "dimensions.h"
1028 #include "paramet.h"
1029 
1030  INTEGER, INTENT(IN) :: timestep
1031 
1032  LOGICAL, SAVE :: first=.true.
1033 ! Identification fichiers et variables NetCDF:
1034  INTEGER, SAVE :: ncidu,varidu,ncidv,varidv,ncidQ
1035  INTEGER, SAVE :: varidQ,ncidt,varidt,ncidps,varidps
1036  INTEGER :: ncidpl,varidpl,varidap,varidbp
1037 ! Variables auxiliaires NetCDF:
1038  INTEGER, DIMENSION(4) :: start,count
1039  INTEGER :: status,rcode
1040 
1041  CHARACTER (len = 80) :: abort_message
1042  CHARACTER (len = 20) :: modname = 'guide_read'
1043 ! -----------------------------------------------------------------
1044 ! Premier appel: initialisation de la lecture des fichiers
1045 ! -----------------------------------------------------------------
1046  if (first) then
1047  ncidpl=-99
1048  print*,'Guide: ouverture des fichiers guidage '
1049 ! Niveaux de pression si non constants
1050  if (guide_modele) then
1051  print *,'Lecture du guidage sur niveaux modele'
1052  rcode = nf90_open('apbp.nc', nf90_nowrite, ncidpl)
1053  IF (rcode.NE.nf_noerr) THEN
1054  print *,'Guide: probleme -> pas de fichier apbp.nc'
1055  CALL abort_gcm(modname,abort_message,1)
1056  ENDIF
1057  rcode = nf90_inq_varid(ncidpl, 'AP', varidap)
1058  IF (rcode.NE.nf_noerr) THEN
1059  print *,'Guide: probleme -> pas de variable AP, fichier apbp.nc'
1060  CALL abort_gcm(modname,abort_message,1)
1061  ENDIF
1062  rcode = nf90_inq_varid(ncidpl, 'BP', varidbp)
1063  IF (rcode.NE.nf_noerr) THEN
1064  print *,'Guide: probleme -> pas de variable BP, fichier apbp.nc'
1065  CALL abort_gcm(modname,abort_message,1)
1066  ENDIF
1067  print*,'ncidpl,varidap',ncidpl,varidap
1068  endif
1069 ! Vent zonal
1070  if (guide_u) then
1071  rcode = nf90_open('u.nc', nf90_nowrite, ncidu)
1072  IF (rcode.NE.nf_noerr) THEN
1073  print *,'Guide: probleme -> pas de fichier u.nc'
1074  CALL abort_gcm(modname,abort_message,1)
1075  ENDIF
1076  rcode = nf90_inq_varid(ncidu, 'UWND', varidu)
1077  IF (rcode.NE.nf_noerr) THEN
1078  print *,'Guide: probleme -> pas de variable UWND, fichier u.nc'
1079  CALL abort_gcm(modname,abort_message,1)
1080  ENDIF
1081  print*,'ncidu,varidu',ncidu,varidu
1082  if (ncidpl.eq.-99) ncidpl=ncidu
1083  endif
1084 ! Vent meridien
1085  if (guide_v) then
1086  rcode = nf90_open('v.nc', nf90_nowrite, ncidv)
1087  IF (rcode.NE.nf_noerr) THEN
1088  print *,'Guide: probleme -> pas de fichier v.nc'
1089  CALL abort_gcm(modname,abort_message,1)
1090  ENDIF
1091  rcode = nf90_inq_varid(ncidv, 'VWND', varidv)
1092  IF (rcode.NE.nf_noerr) THEN
1093  print *,'Guide: probleme -> pas de variable VWND, fichier v.nc'
1094  CALL abort_gcm(modname,abort_message,1)
1095  ENDIF
1096  print*,'ncidv,varidv',ncidv,varidv
1097  if (ncidpl.eq.-99) ncidpl=ncidv
1098  endif
1099 ! Temperature
1100  if (guide_t) then
1101  rcode = nf90_open('T.nc', nf90_nowrite, ncidt)
1102  IF (rcode.NE.nf_noerr) THEN
1103  print *,'Guide: probleme -> pas de fichier T.nc'
1104  CALL abort_gcm(modname,abort_message,1)
1105  ENDIF
1106  rcode = nf90_inq_varid(ncidt, 'AIR', varidt)
1107  IF (rcode.NE.nf_noerr) THEN
1108  print *,'Guide: probleme -> pas de variable AIR, fichier T.nc'
1109  CALL abort_gcm(modname,abort_message,1)
1110  ENDIF
1111  print*,'ncidT,varidT',ncidt,varidt
1112  if (ncidpl.eq.-99) ncidpl=ncidt
1113  endif
1114 ! Humidite
1115  if (guide_q) then
1116  rcode = nf90_open('hur.nc', nf90_nowrite, ncidq)
1117  IF (rcode.NE.nf_noerr) THEN
1118  print *,'Guide: probleme -> pas de fichier hur.nc'
1119  CALL abort_gcm(modname,abort_message,1)
1120  ENDIF
1121  rcode = nf90_inq_varid(ncidq, 'RH', varidq)
1122  IF (rcode.NE.nf_noerr) THEN
1123  print *,'Guide: probleme -> pas de variable RH, fichier hur.nc'
1124  CALL abort_gcm(modname,abort_message,1)
1125  ENDIF
1126  print*,'ncidQ,varidQ',ncidq,varidq
1127  if (ncidpl.eq.-99) ncidpl=ncidq
1128  endif
1129 ! Pression de surface
1130  if ((guide_p).OR.(guide_modele)) then
1131  rcode = nf90_open('ps.nc', nf90_nowrite, ncidps)
1132  IF (rcode.NE.nf_noerr) THEN
1133  print *,'Guide: probleme -> pas de fichier ps.nc'
1134  CALL abort_gcm(modname,abort_message,1)
1135  ENDIF
1136  rcode = nf90_inq_varid(ncidps, 'SP', varidps)
1137  IF (rcode.NE.nf_noerr) THEN
1138  print *,'Guide: probleme -> pas de variable SP, fichier ps.nc'
1139  CALL abort_gcm(modname,abort_message,1)
1140  ENDIF
1141  print*,'ncidps,varidps',ncidps,varidps
1142  endif
1143 ! Coordonnee verticale
1144  if (.not.guide_modele) then
1145  rcode = nf90_inq_varid(ncidpl, 'LEVEL', varidpl)
1146  IF (rcode.NE.0) rcode = nf90_inq_varid(ncidpl, 'PRESSURE', varidpl)
1147  print*,'ncidpl,varidpl',ncidpl,varidpl
1148  endif
1149 ! Coefs ap, bp pour calcul de la pression aux differents niveaux
1150  if (guide_modele) then
1151 #ifdef NC_DOUBLE
1152  status=nf_get_vara_double(ncidpl,varidap,1,nlevnc,apnc)
1153  status=nf_get_vara_double(ncidpl,varidbp,1,nlevnc,bpnc)
1154 #else
1155  status=nf_get_vara_real(ncidpl,varidap,1,nlevnc,apnc)
1156  status=nf_get_vara_real(ncidpl,varidbp,1,nlevnc,bpnc)
1157 #endif
1158  else
1159 #ifdef NC_DOUBLE
1160  status=nf_get_vara_double(ncidpl,varidpl,1,nlevnc,apnc)
1161 #else
1162  status=nf_get_vara_real(ncidpl,varidpl,1,nlevnc,apnc)
1163 #endif
1164  apnc=apnc*100.! conversion en Pascals
1165  bpnc(:)=0.
1166  endif
1167  first=.false.
1168  endif ! (first)
1169 
1170 ! -----------------------------------------------------------------
1171 ! lecture des champs u, v, T, Q, ps
1172 ! -----------------------------------------------------------------
1173 
1174 ! dimensions pour les champs scalaires et le vent zonal
1175  start(1)=1
1176  start(2)=1
1177  start(3)=1
1178  start(4)=timestep
1179 
1180  count(1)=iip1
1181  count(2)=jjp1
1182  count(3)=nlevnc
1183  count(4)=1
1184 
1185 ! Vent zonal
1186  if (guide_u) then
1187 #ifdef NC_DOUBLE
1188  status=nf_get_vara_double(ncidu,varidu,start,count,unat2)
1189 #else
1190  status=nf_get_vara_real(ncidu,varidu,start,count,unat2)
1191 #endif
1192  IF (invert_y) THEN
1193  CALL invert_lat(iip1,jjp1,nlevnc,unat2)
1194  ENDIF
1195  endif
1196 
1197 ! Temperature
1198  if (guide_t) then
1199 #ifdef NC_DOUBLE
1200  status=nf_get_vara_double(ncidt,varidt,start,count,tnat2)
1201 #else
1202  status=nf_get_vara_real(ncidt,varidt,start,count,tnat2)
1203 #endif
1204  IF (invert_y) THEN
1205  CALL invert_lat(iip1,jjp1,nlevnc,tnat2)
1206  ENDIF
1207  endif
1208 
1209 ! Humidite
1210  if (guide_q) then
1211 #ifdef NC_DOUBLE
1212  status=nf_get_vara_double(ncidq,varidq,start,count,qnat2)
1213 #else
1214  status=nf_get_vara_real(ncidq,varidq,start,count,qnat2)
1215 #endif
1216  IF (invert_y) THEN
1217  CALL invert_lat(iip1,jjp1,nlevnc,qnat2)
1218  ENDIF
1219 
1220  endif
1221 
1222 ! Vent meridien
1223  if (guide_v) then
1224  count(2)=jjm
1225 #ifdef NC_DOUBLE
1226  status=nf_get_vara_double(ncidv,varidv,start,count,vnat2)
1227 #else
1228  status=nf_get_vara_real(ncidv,varidv,start,count,vnat2)
1229 #endif
1230  IF (invert_y) THEN
1231  CALL invert_lat(iip1,jjm,nlevnc,vnat2)
1232  ENDIF
1233  endif
1234 
1235 ! Pression de surface
1236  if ((guide_p).OR.(guide_modele)) then
1237  start(3)=timestep
1238  start(4)=0
1239  count(2)=jjp1
1240  count(3)=1
1241  count(4)=0
1242 #ifdef NC_DOUBLE
1243  status=nf_get_vara_double(ncidps,varidps,start,count,psnat2)
1244 #else
1245  status=nf_get_vara_real(ncidps,varidps,start,count,psnat2)
1246 #endif
1247  IF (invert_y) THEN
1248  CALL invert_lat(iip1,jjp1,1,psnat2)
1249  ENDIF
1250  endif
1251 
1252  END SUBROUTINE guide_read
1253 
1254 !=======================================================================
1255  SUBROUTINE guide_read2d(timestep)
1256 
1257  IMPLICIT NONE
1258 
1259 #include "netcdf.inc"
1260 #include "dimensions.h"
1261 #include "paramet.h"
1262 
1263  INTEGER, INTENT(IN) :: timestep
1264 
1265  LOGICAL, SAVE :: first=.true.
1266 ! Identification fichiers et variables NetCDF:
1267  INTEGER, SAVE :: ncidu,varidu,ncidv,varidv,ncidQ
1268  INTEGER, SAVE :: varidQ,ncidt,varidt,ncidps,varidps
1269  INTEGER :: ncidpl,varidpl,varidap,varidbp
1270 ! Variables auxiliaires NetCDF:
1271  INTEGER, DIMENSION(4) :: start,count
1272  INTEGER :: status,rcode
1273 ! Variables for 3D extension:
1274  REAL, DIMENSION (jjp1,llm) :: zu
1275  REAL, DIMENSION (jjm,llm) :: zv
1276  INTEGER :: i
1277 
1278  CHARACTER (len = 80) :: abort_message
1279  CHARACTER (len = 20) :: modname = 'guide_read2D'
1280 ! -----------------------------------------------------------------
1281 ! Premier appel: initialisation de la lecture des fichiers
1282 ! -----------------------------------------------------------------
1283  if (first) then
1284  ncidpl=-99
1285  print*,'Guide: ouverture des fichiers guidage '
1286 ! Niveaux de pression si non constants
1287  if (guide_modele) then
1288  print *,'Lecture du guidage sur niveaux modele'
1289  rcode = nf90_open('apbp.nc', nf90_nowrite, ncidpl)
1290  IF (rcode.NE.nf_noerr) THEN
1291  print *,'Guide: probleme -> pas de fichier apbp.nc'
1292  CALL abort_gcm(modname,abort_message,1)
1293  ENDIF
1294  rcode = nf90_inq_varid(ncidpl, 'AP', varidap)
1295  IF (rcode.NE.nf_noerr) THEN
1296  print *,'Guide: probleme -> pas de variable AP, fichier apbp.nc'
1297  CALL abort_gcm(modname,abort_message,1)
1298  ENDIF
1299  rcode = nf90_inq_varid(ncidpl, 'BP', varidbp)
1300  IF (rcode.NE.nf_noerr) THEN
1301  print *,'Guide: probleme -> pas de variable BP, fichier apbp.nc'
1302  CALL abort_gcm(modname,abort_message,1)
1303  ENDIF
1304  print*,'ncidpl,varidap',ncidpl,varidap
1305  endif
1306 ! Vent zonal
1307  if (guide_u) then
1308  rcode = nf90_open('u.nc', nf90_nowrite, ncidu)
1309  IF (rcode.NE.nf_noerr) THEN
1310  print *,'Guide: probleme -> pas de fichier u.nc'
1311  CALL abort_gcm(modname,abort_message,1)
1312  ENDIF
1313  rcode = nf90_inq_varid(ncidu, 'UWND', varidu)
1314  IF (rcode.NE.nf_noerr) THEN
1315  print *,'Guide: probleme -> pas de variable UWND, fichier u.nc'
1316  CALL abort_gcm(modname,abort_message,1)
1317  ENDIF
1318  print*,'ncidu,varidu',ncidu,varidu
1319  if (ncidpl.eq.-99) ncidpl=ncidu
1320  endif
1321 ! Vent meridien
1322  if (guide_v) then
1323  rcode = nf90_open('v.nc', nf90_nowrite, ncidv)
1324  IF (rcode.NE.nf_noerr) THEN
1325  print *,'Guide: probleme -> pas de fichier v.nc'
1326  CALL abort_gcm(modname,abort_message,1)
1327  ENDIF
1328  rcode = nf90_inq_varid(ncidv, 'VWND', varidv)
1329  IF (rcode.NE.nf_noerr) THEN
1330  print *,'Guide: probleme -> pas de variable VWND, fichier v.nc'
1331  CALL abort_gcm(modname,abort_message,1)
1332  ENDIF
1333  print*,'ncidv,varidv',ncidv,varidv
1334  if (ncidpl.eq.-99) ncidpl=ncidv
1335  endif
1336 ! Temperature
1337  if (guide_t) then
1338  rcode = nf90_open('T.nc', nf90_nowrite, ncidt)
1339  IF (rcode.NE.nf_noerr) THEN
1340  print *,'Guide: probleme -> pas de fichier T.nc'
1341  CALL abort_gcm(modname,abort_message,1)
1342  ENDIF
1343  rcode = nf90_inq_varid(ncidt, 'AIR', varidt)
1344  IF (rcode.NE.nf_noerr) THEN
1345  print *,'Guide: probleme -> pas de variable AIR, fichier T.nc'
1346  CALL abort_gcm(modname,abort_message,1)
1347  ENDIF
1348  print*,'ncidT,varidT',ncidt,varidt
1349  if (ncidpl.eq.-99) ncidpl=ncidt
1350  endif
1351 ! Humidite
1352  if (guide_q) then
1353  rcode = nf90_open('hur.nc', nf90_nowrite, ncidq)
1354  IF (rcode.NE.nf_noerr) THEN
1355  print *,'Guide: probleme -> pas de fichier hur.nc'
1356  CALL abort_gcm(modname,abort_message,1)
1357  ENDIF
1358  rcode = nf90_inq_varid(ncidq, 'RH', varidq)
1359  IF (rcode.NE.nf_noerr) THEN
1360  print *,'Guide: probleme -> pas de variable RH, fichier hur.nc'
1361  CALL abort_gcm(modname,abort_message,1)
1362  ENDIF
1363  print*,'ncidQ,varidQ',ncidq,varidq
1364  if (ncidpl.eq.-99) ncidpl=ncidq
1365  endif
1366 ! Pression de surface
1367  if ((guide_p).OR.(guide_modele)) then
1368  rcode = nf90_open('ps.nc', nf90_nowrite, ncidps)
1369  IF (rcode.NE.nf_noerr) THEN
1370  print *,'Guide: probleme -> pas de fichier ps.nc'
1371  CALL abort_gcm(modname,abort_message,1)
1372  ENDIF
1373  rcode = nf90_inq_varid(ncidps, 'SP', varidps)
1374  IF (rcode.NE.nf_noerr) THEN
1375  print *,'Guide: probleme -> pas de variable SP, fichier ps.nc'
1376  CALL abort_gcm(modname,abort_message,1)
1377  ENDIF
1378  print*,'ncidps,varidps',ncidps,varidps
1379  endif
1380 ! Coordonnee verticale
1381  if (.not.guide_modele) then
1382  rcode = nf90_inq_varid(ncidpl, 'LEVEL', varidpl)
1383  IF (rcode.NE.0) rcode = nf90_inq_varid(ncidpl, 'PRESSURE', varidpl)
1384  print*,'ncidpl,varidpl',ncidpl,varidpl
1385  endif
1386 ! Coefs ap, bp pour calcul de la pression aux differents niveaux
1387  if (guide_modele) then
1388 #ifdef NC_DOUBLE
1389  status=nf_get_vara_double(ncidpl,varidap,1,nlevnc,apnc)
1390  status=nf_get_vara_double(ncidpl,varidbp,1,nlevnc,bpnc)
1391 #else
1392  status=nf_get_vara_real(ncidpl,varidap,1,nlevnc,apnc)
1393  status=nf_get_vara_real(ncidpl,varidbp,1,nlevnc,bpnc)
1394 #endif
1395  else
1396 #ifdef NC_DOUBLE
1397  status=nf_get_vara_double(ncidpl,varidpl,1,nlevnc,apnc)
1398 #else
1399  status=nf_get_vara_real(ncidpl,varidpl,1,nlevnc,apnc)
1400 #endif
1401  apnc=apnc*100.! conversion en Pascals
1402  bpnc(:)=0.
1403  endif
1404  first=.false.
1405  endif ! (first)
1406 
1407 ! -----------------------------------------------------------------
1408 ! lecture des champs u, v, T, Q, ps
1409 ! -----------------------------------------------------------------
1410 
1411 ! dimensions pour les champs scalaires et le vent zonal
1412  start(1)=1
1413  start(2)=1
1414  start(3)=1
1415  start(4)=timestep
1416 
1417  count(1)=1
1418  count(2)=jjp1
1419  count(3)=nlevnc
1420  count(4)=1
1421 
1422 ! Vent zonal
1423  if (guide_u) then
1424 #ifdef NC_DOUBLE
1425  status=nf_get_vara_double(ncidu,varidu,start,count,zu)
1426 #else
1427  status=nf_get_vara_real(ncidu,varidu,start,count,zu)
1428 #endif
1429  DO i=1,iip1
1430  unat2(i,:,:)=zu(:,:)
1431  ENDDO
1432 
1433  IF (invert_y) THEN
1434  CALL invert_lat(iip1,jjp1,nlevnc,unat2)
1435  ENDIF
1436 
1437  endif
1438 
1439 ! Temperature
1440  if (guide_t) then
1441 #ifdef NC_DOUBLE
1442  status=nf_get_vara_double(ncidt,varidt,start,count,zu)
1443 #else
1444  status=nf_get_vara_real(ncidt,varidt,start,count,zu)
1445 #endif
1446  DO i=1,iip1
1447  tnat2(i,:,:)=zu(:,:)
1448  ENDDO
1449 
1450  IF (invert_y) THEN
1451  CALL invert_lat(iip1,jjp1,nlevnc,tnat2)
1452  ENDIF
1453 
1454  endif
1455 
1456 ! Humidite
1457  if (guide_q) then
1458 #ifdef NC_DOUBLE
1459  status=nf_get_vara_double(ncidq,varidq,start,count,zu)
1460 #else
1461  status=nf_get_vara_real(ncidq,varidq,start,count,zu)
1462 #endif
1463  DO i=1,iip1
1464  qnat2(i,:,:)=zu(:,:)
1465  ENDDO
1466 
1467  IF (invert_y) THEN
1468  CALL invert_lat(iip1,jjp1,nlevnc,qnat2)
1469  ENDIF
1470 
1471  endif
1472 
1473 ! Vent meridien
1474  if (guide_v) then
1475  count(2)=jjm
1476 #ifdef NC_DOUBLE
1477  status=nf_get_vara_double(ncidv,varidv,start,count,zv)
1478 #else
1479  status=nf_get_vara_real(ncidv,varidv,start,count,zv)
1480 #endif
1481  DO i=1,iip1
1482  vnat2(i,:,:)=zv(:,:)
1483  ENDDO
1484 
1485  IF (invert_y) THEN
1486  CALL invert_lat(iip1,jjm,nlevnc,vnat2)
1487  ENDIF
1488 
1489  endif
1490 
1491 ! Pression de surface
1492  if ((guide_p).OR.(guide_modele)) then
1493  start(3)=timestep
1494  start(4)=0
1495  count(2)=jjp1
1496  count(3)=1
1497  count(4)=0
1498 #ifdef NC_DOUBLE
1499  status=nf_get_vara_double(ncidps,varidps,start,count,zu(:,1))
1500 #else
1501  status=nf_get_vara_real(ncidps,varidps,start,count,zu(:,1))
1502 #endif
1503  DO i=1,iip1
1504  psnat2(i,:)=zu(:,1)
1505  ENDDO
1506 
1507  IF (invert_y) THEN
1508  CALL invert_lat(iip1,jjp1,1,psnat2)
1509  ENDIF
1510 
1511  endif
1512 
1513  END SUBROUTINE guide_read2d
1514 
1515 !=======================================================================
1516  SUBROUTINE guide_out(varname,hsize,vsize,field)
1517 
1518  IMPLICIT NONE
1519 
1520  include "dimensions.h"
1521  include "paramet.h"
1522  include "netcdf.inc"
1523  include "comgeom2.h"
1524  include "comconst.h"
1525  include "comvert.h"
1526 
1527  ! Variables entree
1528  CHARACTER*(*), INTENT(IN) :: varname
1529  INTEGER, INTENT (IN) :: hsize,vsize
1530  REAL, DIMENSION (iip1,hsize,vsize), INTENT(IN) :: field
1531 
1532  ! Variables locales
1533  INTEGER, SAVE :: timestep=0
1534  ! Identites fichier netcdf
1535  INTEGER :: nid, id_lonu, id_lonv, id_latu, id_latv, id_tim, id_lev
1536  INTEGER :: vid_lonu,vid_lonv,vid_latu,vid_latv,vid_cu,vid_cv,vid_lev
1537  INTEGER :: vid_au,vid_av
1538  INTEGER, DIMENSION (3) :: dim3
1539  INTEGER, DIMENSION (4) :: dim4,count,start
1540  INTEGER :: ierr, varid,l
1541  REAL, DIMENSION (iip1,hsize,vsize) :: field2
1542 
1543  print *,'Guide: output timestep',timestep,'var ',varname
1544  IF (timestep.EQ.0) THEN
1545 ! ----------------------------------------------
1546 ! initialisation fichier de sortie
1547 ! ----------------------------------------------
1548 ! Ouverture du fichier
1549  ierr=nf_create("guide_ins.nc",nf_clobber,nid)
1550 ! Definition des dimensions
1551  ierr=nf_def_dim(nid,"LONU",iip1,id_lonu)
1552  ierr=nf_def_dim(nid,"LONV",iip1,id_lonv)
1553  ierr=nf_def_dim(nid,"LATU",jjp1,id_latu)
1554  ierr=nf_def_dim(nid,"LATV",jjm,id_latv)
1555  ierr=nf_def_dim(nid,"LEVEL",llm,id_lev)
1556  ierr=nf_def_dim(nid,"TIME",nf_unlimited,id_tim)
1557 
1558 ! Creation des variables dimensions
1559  ierr=nf_def_var(nid,"LONU",nf_float,1,id_lonu,vid_lonu)
1560  ierr=nf_def_var(nid,"LONV",nf_float,1,id_lonv,vid_lonv)
1561  ierr=nf_def_var(nid,"LATU",nf_float,1,id_latu,vid_latu)
1562  ierr=nf_def_var(nid,"LATV",nf_float,1,id_latv,vid_latv)
1563  ierr=nf_def_var(nid,"LEVEL",nf_float,1,id_lev,vid_lev)
1564  ierr=nf_def_var(nid,"cu",nf_float,2,(/id_lonu,id_latu/),vid_cu)
1565  ierr=nf_def_var(nid,"au",nf_float,2,(/id_lonu,id_latu/),vid_au)
1566  ierr=nf_def_var(nid,"cv",nf_float,2,(/id_lonv,id_latv/),vid_cv)
1567  ierr=nf_def_var(nid,"av",nf_float,2,(/id_lonv,id_latv/),vid_av)
1568 
1569  ierr=nf_enddef(nid)
1570 
1571 ! Enregistrement des variables dimensions
1572 #ifdef NC_DOUBLE
1573  ierr = nf_put_var_double(nid,vid_lonu,rlonu*180./pi)
1574  ierr = nf_put_var_double(nid,vid_lonv,rlonv*180./pi)
1575  ierr = nf_put_var_double(nid,vid_latu,rlatu*180./pi)
1576  ierr = nf_put_var_double(nid,vid_latv,rlatv*180./pi)
1577  ierr = nf_put_var_double(nid,vid_lev,presnivs)
1578  ierr = nf_put_var_double(nid,vid_cu,cu)
1579  ierr = nf_put_var_double(nid,vid_cv,cv)
1580  ierr = nf_put_var_double(nid,vid_au,alpha_u)
1581  ierr = nf_put_var_double(nid,vid_av,alpha_v)
1582 #else
1583  ierr = nf_put_var_real(nid,vid_lonu,rlonu*180./pi)
1584  ierr = nf_put_var_real(nid,vid_lonv,rlonv*180./pi)
1585  ierr = nf_put_var_real(nid,vid_latu,rlatu*180./pi)
1586  ierr = nf_put_var_real(nid,vid_latv,rlatv*180./pi)
1587  ierr = nf_put_var_real(nid,vid_lev,presnivs)
1588  ierr = nf_put_var_real(nid,vid_cu,cu)
1589  ierr = nf_put_var_real(nid,vid_cv,cv)
1590  ierr = nf_put_var_real(nid,vid_au,alpha_u)
1591  ierr = nf_put_var_real(nid,vid_av,alpha_v)
1592 #endif
1593 ! --------------------------------------------------------------------
1594 ! Cr�ation des variables sauvegard�es
1595 ! --------------------------------------------------------------------
1596  ierr = nf_redef(nid)
1597 ! Surface pressure (GCM)
1598  dim3=(/id_lonv,id_latu,id_tim/)
1599  ierr = nf_def_var(nid,"SP",nf_float,3,dim3,varid)
1600 ! Surface pressure (guidage)
1601  IF (guide_p) THEN
1602  dim3=(/id_lonv,id_latu,id_tim/)
1603  ierr = nf_def_var(nid,"ps",nf_float,3,dim3,varid)
1604  ENDIF
1605 ! Zonal wind
1606  IF (guide_u) THEN
1607  dim4=(/id_lonu,id_latu,id_lev,id_tim/)
1608  ierr = nf_def_var(nid,"u",nf_float,4,dim4,varid)
1609  ierr = nf_def_var(nid,"ua",nf_float,4,dim4,varid)
1610  ierr = nf_def_var(nid,"ucov",nf_float,4,dim4,varid)
1611  ENDIF
1612 ! Merid. wind
1613  IF (guide_v) THEN
1614  dim4=(/id_lonv,id_latv,id_lev,id_tim/)
1615  ierr = nf_def_var(nid,"v",nf_float,4,dim4,varid)
1616  ierr = nf_def_var(nid,"va",nf_float,4,dim4,varid)
1617  ierr = nf_def_var(nid,"vcov",nf_float,4,dim4,varid)
1618  ENDIF
1619 ! Pot. Temperature
1620  IF (guide_t) THEN
1621  dim4=(/id_lonv,id_latu,id_lev,id_tim/)
1622  ierr = nf_def_var(nid,"teta",nf_float,4,dim4,varid)
1623  ENDIF
1624 ! Specific Humidity
1625  IF (guide_q) THEN
1626  dim4=(/id_lonv,id_latu,id_lev,id_tim/)
1627  ierr = nf_def_var(nid,"q",nf_float,4,dim4,varid)
1628  ENDIF
1629 
1630  ierr = nf_enddef(nid)
1631  ierr = nf_close(nid)
1632  ENDIF ! timestep=0
1633 
1634 ! --------------------------------------------------------------------
1635 ! Enregistrement du champ
1636 ! --------------------------------------------------------------------
1637  ierr=nf_open("guide_ins.nc",nf_write,nid)
1638 
1639  IF (varname=="SP") timestep=timestep+1
1640 
1641  ierr = nf_inq_varid(nid,varname,varid)
1642  SELECT CASE (varname)
1643  CASE ("SP","ps")
1644  start=(/1,1,timestep,0/)
1645  count=(/iip1,jjp1,1,0/)
1646  CASE ("v","va","vcov")
1647  start=(/1,1,1,timestep/)
1648  count=(/iip1,jjm,llm,1/)
1649  CASE DEFAULT
1650  start=(/1,1,1,timestep/)
1651  count=(/iip1,jjp1,llm,1/)
1652  END SELECT
1653 
1654  SELECT CASE (varname)
1655  CASE("u","ua")
1656  DO l=1,llm ; field2(:,2:jjm,l)=field(:,2:jjm,l)/cu(:,2:jjm) ; ENDDO
1657  field2(:,1,:)=0. ; field2(:,jjp1,:)=0.
1658  CASE("v","va")
1659  DO l=1,llm ; field2(:,:,l)=field(:,:,l)/cv(:,:) ; ENDDO
1660  CASE DEFAULT
1661  field2=field
1662  END SELECT
1663 
1664 
1665 #ifdef NC_DOUBLE
1666  ierr = nf_put_vara_double(nid,varid,start,count,field2)
1667 #else
1668  ierr = nf_put_vara_real(nid,varid,start,count,field2)
1669 #endif
1670 
1671  ierr = nf_close(nid)
1672 
1673  END SUBROUTINE guide_out
1674 
1675 
1676 !===========================================================================
1677  subroutine correctbid(iim,nl,x)
1678  integer iim,nl
1679  real x(iim+1,nl)
1680  integer i,l
1681  real zz
1682 
1683  do l=1,nl
1684  do i=2,iim-1
1685  if(abs(x(i,l)).gt.1.e10) then
1686  zz=0.5*(x(i-1,l)+x(i+1,l))
1687  print*,'correction ',i,l,x(i,l),zz
1688  x(i,l)=zz
1689  endif
1690  enddo
1691  enddo
1692  return
1693  end subroutine correctbid
1694 
1695 !===========================================================================
1696 END MODULE guide_mod
guide_mod::unat2
real, dimension(:,:,:), allocatable, save, private unat2
Definition: guide_mod.F90:47
ip1jmp1
!$Header llmm1 INTEGER ip1jmp1
Definition: paramet.h:14
guide_mod::unat1
real, dimension(:,:,:), allocatable, save, private unat1
Definition: guide_mod.F90:47
guide_mod::lat_max_g
real, save, private lat_max_g
Definition: guide_mod.F90:35
grossismx
!$Header!c!c!c include serre h!c REAL && grossismx
Definition: serre.h:8
exner_milieu_m
Definition: exner_milieu_m.F90:1
write_field::writefield
Definition: write_field.F90:15
guide_mod::guide_addfield
subroutine guide_addfield(hsize, vsize, field, alpha)
Definition: guide_mod.F90:522
guide_mod::tgui2
real, dimension(:,:), allocatable, save, private tgui2
Definition: guide_mod.F90:56
guide_mod::vnat2
real, dimension(:,:,:), allocatable, save, private vnat2
Definition: guide_mod.F90:48
getparam::fin_getparam
subroutine fin_getparam
Definition: getparam.F90:29
pres2lev_mod::pres2lev
subroutine pres2lev(varo, varn, lmo, lmn, po, pn, ni, nj, ok_invertp)
Definition: pres2lev_mod.F90:9
preff
!$Id preff
Definition: comvert.h:8
guide_mod::guide_main
subroutine guide_main(itau, ucov, vcov, teta, q, masse, ps)
Definition: guide_mod.F90:315
guide_mod::lat_min_g
real, save, private lat_min_g
Definition: guide_mod.F90:35
bp
!$Id bp(llm+1)
kappa
!$Id mode_top_bound COMMON comconstr kappa
Definition: comconst.h:7
guide_mod::ugui1
real, dimension(:,:), allocatable, save, private ugui1
Definition: guide_mod.F90:54
guide_mod::tgui1
real, dimension(:,:), allocatable, save, private tgui1
Definition: guide_mod.F90:56
exner_hyb_m::exner_hyb
subroutine exner_hyb(ngrid, ps, p, pks, pk, pkf)
Definition: exner_hyb_m.F90:8
guide_mod::tau_lat
real, save, private tau_lat
Definition: guide_mod.F90:37
guide_mod::ini_anal
logical, save, private ini_anal
Definition: guide_mod.F90:26
clon
!$Header!c!c!c include serre h!c REAL clon
Definition: serre.h:8
invert_lat
subroutine invert_lat(xsize, ysize, vsize, field)
Definition: invert_lat.F90:3
guide_mod::tnat1
real, dimension(:,:,:), allocatable, save, private tnat1
Definition: guide_mod.F90:49
exner_milieu_m::exner_milieu
subroutine exner_milieu(ngrid, ps, p, pks, pk, pkf)
Definition: exner_milieu_m.F90:8
abort_gcm
subroutine abort_gcm(modname, message, ierr)
Definition: abort_gcm.F:7
llm
!$Id Turb_fcg_gcssold get_uvd hqturb_gcssold endif!large scale llm day day1 day day1 *dt_toga endif!time annee_ref dt_toga u_toga vq_toga w_prof vq_prof llm day day1 day day1 *dt_dice endif!time annee_ref dt_dice swup_dice vg_dice omega_dice tg_prof vg_profd w_profd omega_profd!do llm!print llm l llm
Definition: 1D_interp_cases.h:103
guide_mod::tau_min_v
real, save, private tau_min_v
Definition: guide_mod.F90:30
guide_mod::invert_p
logical, save, private invert_p
Definition: guide_mod.F90:26
pi
!$Id mode_top_bound COMMON comconstr && pi
Definition: comconst.h:7
aire
!$Header!CDK comgeom COMMON comgeom aire
Definition: comgeom.h:25
getparam::ini_getparam
subroutine ini_getparam(fichier)
Definition: getparam.F90:21
rlatu
!$Header!CDK comgeom COMMON comgeom rlatu
Definition: comgeom.h:25
control_mod::day_step
integer, save day_step
Definition: control_mod.F90:15
guide_mod::invert_y
logical, save, private invert_y
Definition: guide_mod.F90:26
presnivs
!$Id presnivs(llm)
guide_mod::guide_hr
logical, save, private guide_hr
Definition: guide_mod.F90:24
guide_mod::tau_lon
real, save, private tau_lon
Definition: guide_mod.F90:37
guide_mod::lon_max_g
real, save, private lon_max_g
Definition: guide_mod.F90:36
guide_mod::alpha_q
real, dimension(:), allocatable, save, private alpha_q
Definition: guide_mod.F90:40
guide_mod::tau_max_t
real, save, private tau_max_t
Definition: guide_mod.F90:31
guide_mod::tau2alpha
subroutine tau2alpha(typ, pim, pjm, factt, taumin, taumax, alpha)
Definition: guide_mod.F90:854
guide_mod::guide_read2d
subroutine guide_read2d(timestep)
Definition: guide_mod.F90:1256
ip1jm
!$Header llmm1 INTEGER ip1jm
Definition: paramet.h:14
pression
subroutine pression(ngrid, ap, bp, ps, p)
Definition: pression.F90:2
massbar
subroutine massbar(masse, massebx, masseby)
Definition: massbar.F90:2
false
!$Id itapm1 ENDIF!IM on interpole les champs sur les niveaux STD de pression!IM a chaque pas de temps de la physique c!positionnement de l argument logique a false c!pour ne pas recalculer deux fois la meme chose!c!a cet effet un appel a plevel_new a ete deplace c!a la fin de la serie d appels c!la boucle DO nlevSTD a ete internalisee c!dans d ou la creation de cette routine c c!CALL false
Definition: calcul_STDlev.h:26
guide_mod::tau_max_p
real, save, private tau_max_p
Definition: guide_mod.F90:33
guide_mod::alpha_v
real, dimension(:), allocatable, save, private alpha_v
Definition: guide_mod.F90:39
q_sat
subroutine q_sat(np, temp, pres, qsat)
Definition: q_sat.F:8
guide_mod::alpha_pcor
real, dimension(:), allocatable, save, private alpha_pcor
Definition: guide_mod.F90:41
guide_mod::tau_max_v
real, save, private tau_max_v
Definition: guide_mod.F90:30
aireu
!$Header!CDK comgeom COMMON comgeom aireu
Definition: comgeom.h:25
guide_mod::guide_sav
logical, save, private guide_sav
Definition: guide_mod.F90:27
jjp1
!$Header jjp1
Definition: paramet.h:14
guide_mod::guide_u
logical, save, private guide_u
Definition: guide_mod.F90:23
guide_mod::guide_bl
logical, save, private guide_bl
Definition: guide_mod.F90:25
guide_mod::bpnc
real, dimension(:), allocatable, save, private bpnc
Definition: guide_mod.F90:52
cpp
!$Id mode_top_bound COMMON comconstr cpp
Definition: comconst.h:7
guide_mod::psnat2
real, dimension(:,:), allocatable, save, private psnat2
Definition: guide_mod.F90:51
guide_mod::iguide_read
integer, save, private iguide_read
Definition: guide_mod.F90:21
guide_mod::gamma4
logical, save, private gamma4
Definition: guide_mod.F90:25
guide_mod::tau_min_p
real, save, private tau_min_p
Definition: guide_mod.F90:33
guide_mod::alpha_p
real, dimension(:), allocatable, save, private alpha_p
Definition: guide_mod.F90:41
guide_mod::guide_read
subroutine guide_read(timestep)
Definition: guide_mod.F90:1023
guide_mod::guide_init
subroutine guide_init
Definition: guide_mod.F90:64
guide_mod::guide_2d
logical, save, private guide_2d
Definition: guide_mod.F90:27
daysec
!$Id mode_top_bound COMMON comconstr daysec
Definition: comconst.h:7
guide_mod::psgui2
real, dimension(:), allocatable, save, private psgui2
Definition: guide_mod.F90:58
rlonu
!$Header!CDK comgeom COMMON comgeom rlonu
Definition: comgeom.h:25
rlatv
!$Header!CDK comgeom COMMON comgeom rlatv
Definition: comgeom.h:25
exner_hyb_m
Definition: exner_hyb_m.F90:1
write_field
Definition: write_field.F90:4
guide_mod::alpha_u
real, dimension(:), allocatable, save, private alpha_u
Definition: guide_mod.F90:39
guide_mod::tau_min_q
real, save, private tau_min_q
Definition: guide_mod.F90:32
guide_mod::ugui2
real, dimension(:,:), allocatable, save, private ugui2
Definition: guide_mod.F90:54
guide_mod::guide_v
logical, save, private guide_v
Definition: guide_mod.F90:23
true
!$Id itapm1 ENDIF!IM on interpole les champs sur les niveaux STD de pression!IM a chaque pas de temps de la physique c!positionnement de l argument logique a false c!pour ne pas recalculer deux fois la meme chose!c!a cet effet un appel a plevel_new a ete deplace c!a la fin de la serie d appels c!la boucle DO nlevSTD a ete internalisee c!dans d ou la creation de cette routine c c!CALL ulevSTD CALL &zphi philevSTD CALL &zx_rh rhlevSTD!DO klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon klev DO klon du jour ou toutes les read_climoz CALL true
Definition: calcul_STDlev.h:139
guide_mod::tnat2
real, dimension(:,:,:), allocatable, save, private tnat2
Definition: guide_mod.F90:49
guide_mod
Definition: guide_mod.F90:4
cu
!$Id!Parameters for parameters that control the rate of approach!to quasi equilibrium noff nlm real tlcrit real entp real sigd real coeffs real dtmax real cu real betad real damp real delta COMMON cvparam nlm tlcrit sigd coeffs cu
Definition: cvparam.h:12
coordij
subroutine coordij(lon, lat, ilon, jlat)
Definition: coordij.F:5
guide_mod::tau_min_t
real, save, private tau_min_t
Definition: guide_mod.F90:31
guide_mod::guide_zon
logical, save, private guide_zon
Definition: guide_mod.F90:25
guide_mod::guide_add
logical, save, private guide_add
Definition: guide_mod.F90:25
grossismy
!$Header!c!c!c include serre h!c REAL grossismy
Definition: serre.h:8
guide_mod::guide_interp
subroutine guide_interp(psi, teta)
Definition: guide_mod.F90:606
guide_mod::vgui1
real, dimension(:,:), allocatable, save, private vgui1
Definition: guide_mod.F90:55
guide_mod::tau_max_u
real, save, private tau_max_u
Definition: guide_mod.F90:29
guide_mod::guide_t
logical, save, private guide_t
Definition: guide_mod.F90:23
dtvr
!$Id mode_top_bound COMMON comconstr dtvr
Definition: comconst.h:7
guide_mod::alpha_t
real, dimension(:), allocatable, save, private alpha_t
Definition: guide_mod.F90:40
guide_mod::vnat1
real, dimension(:,:,:), allocatable, save, private vnat1
Definition: guide_mod.F90:48
guide_mod::qnat2
real, dimension(:,:,:), allocatable, save, private qnat2
Definition: guide_mod.F90:50
pres2lev_mod
Definition: pres2lev_mod.F90:3
guide_mod::vgui2
real, dimension(:,:), allocatable, save, private vgui2
Definition: guide_mod.F90:55
iim
c c zjulian c cym CALL iim cym klev iim
Definition: ini_bilKP_ave.h:24
guide_mod::tau_min_u
real, save, private tau_min_u
Definition: guide_mod.F90:29
control_mod::iperiod
integer, save iperiod
Definition: control_mod.F90:16
guide_mod::guide_q
logical, save, private guide_q
Definition: guide_mod.F90:23
guide_mod::lon_min_g
real, save, private lon_min_g
Definition: guide_mod.F90:36
clat
!$Header!c!c!c include serre h!c REAL clat
Definition: serre.h:8
guide_mod::qnat1
real, dimension(:,:,:), allocatable, save, private qnat1
Definition: guide_mod.F90:50
guide_mod::apnc
real, dimension(:), allocatable, save, private apnc
Definition: guide_mod.F90:52
guide_mod::psnat1
real, dimension(:,:), allocatable, save, private psnat1
Definition: guide_mod.F90:51
massdair
subroutine massdair(p, masse)
Definition: massdair.F:5
guide_mod::psgui1
real, dimension(:), allocatable, save, private psgui1
Definition: guide_mod.F90:58
cv
!$Header!CDK comgeom COMMON comgeom cv
Definition: comgeom.h:25
guide_mod::guide_teta
logical, save, private guide_teta
Definition: guide_mod.F90:24
guide_mod::guide_out
subroutine guide_out(varname, hsize, vsize, field)
Definition: guide_mod.F90:1517
guide_mod::guide_modele
logical, save, private guide_modele
Definition: guide_mod.F90:26
guide_mod::qgui1
real, dimension(:,:), allocatable, save, private qgui1
Definition: guide_mod.F90:57
guide_mod::tau_max_q
real, save, private tau_max_q
Definition: guide_mod.F90:32
guide_mod::guide_reg
logical, save, private guide_reg
Definition: guide_mod.F90:25
guide_mod::iguide_int
integer, save, private iguide_int
Definition: guide_mod.F90:21
guide_mod::guide_zonave
subroutine guide_zonave(typ, hsize, vsize, field)
Definition: guide_mod.F90:543
guide_mod::nlevnc
integer, save, private nlevnc
Definition: guide_mod.F90:22
guide_mod::guide_p
logical, save, private guide_p
Definition: guide_mod.F90:23
guide_mod::qgui2
real, dimension(:,:), allocatable, save, private qgui2
Definition: guide_mod.F90:57
control_mod
Definition: control_mod.F90:5
guide_mod::correctbid
subroutine correctbid(iim, nl, x)
Definition: guide_mod.F90:1678
guide_mod::iguide_sav
integer, save, private iguide_sav
Definition: guide_mod.F90:21
airev
!$Header!CDK comgeom COMMON comgeom airev
Definition: comgeom.h:25
rlonv
!$Header!CDK comgeom COMMON comgeom rlonv
Definition: comgeom.h:25