/[lmdze]/trunk/Sources/dyn3d/initdynav.f
ViewVC logotype

Contents of /trunk/Sources/dyn3d/initdynav.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 66 - (show annotations)
Thu Sep 20 13:00:41 2012 UTC (11 years, 8 months ago) by guez
Original Path: trunk/libf/bibio/initdynav.f90
File size: 2965 byte(s)
Changed name of module "comvert" to "disvert_m". Changed constant
1. to 0.3 in vertical sampling "strato".

1 module initdynav_m
2
3 implicit none
4
5 integer histaveid
6
7 contains
8
9 subroutine initdynav(day0, anne0, tstep, nq, t_ops, t_wrt)
10
11 ! From initdynav.F, version 1.1.1.1, 2004/05/19 12:53:05
12 ! L. Fairhead, LMD
13 ! Routine d'initialisation des écritures des fichiers histoires LMDZ
14 ! au format IOIPSL. Initialisation du fichier histoire moyenne.
15
16 use calendar, ONLY: ymds2ju
17 USE disvert_m, ONLY: nivsigs
18 USE comgeom, ONLY: rlatu, rlonv
19 USE dimens_m, ONLY: llm
20 USE histbeg_totreg_m, ONLY: histbeg_totreg
21 USE histdef_m, ONLY: histdef
22 USE histend_m, ONLY: histend
23 USE histvert_m, ONLY: histvert
24 USE iniadvtrac_m, ONLY: ttext
25 USE nr_util, ONLY: pi
26 USE paramet_m, ONLY: iip1, jjp1
27 USE temps, ONLY: itau_dyn
28
29 integer, intent(in):: day0, anne0 ! date de reference
30 real, intent(in):: tstep ! frequence d'ecriture
31 integer, intent(in):: nq ! nombre de traceurs
32 real, intent(in):: t_ops ! frequence de l'operation pour IOIPSL
33 real, intent(in):: t_wrt ! frequence d'ecriture sur le fichier
34
35 ! Variables locales
36 integer horiid, zvertiid
37 real julian
38 integer iq
39 integer ii, jj
40
41 !----------------------------------------------------
42
43 print *, "Call sequence information: initdynav"
44
45 CALL ymds2ju(anne0, 1, day0, 0., julian)
46 call histbeg_totreg('dyn_hist_ave.nc', rlonv * 180. / pi, &
47 rlatu * 180. / pi, 1, iip1, 1, jjp1, itau_dyn, julian, tstep, &
48 horiid, histaveid)
49 call histvert(histaveid, 'sigss', 'Niveaux sigma', 'Pa', llm, nivsigs, &
50 zvertiid)
51
52 call histdef(histaveid, 'u', 'vents u scalaires moyennes', 'm/s', iip1, &
53 jjp1, horiid, llm, 1, llm, zvertiid, 'ave(X)', t_ops, t_wrt)
54 call histdef(histaveid, 'v', 'vents v scalaires moyennes', 'm/s', iip1, &
55 jjp1, horiid, llm, 1, llm, zvertiid, 'ave(X)', t_ops, t_wrt)
56 call histdef(histaveid, 'temp', 'temperature moyennee', 'K', iip1, jjp1, &
57 horiid, llm, 1, llm, zvertiid, 'ave(X)', t_ops, t_wrt)
58 call histdef(histaveid, 'theta', 'temperature potentielle', 'K', iip1, &
59 jjp1, horiid, llm, 1, llm, zvertiid, 'ave(X)', t_ops, t_wrt)
60 call histdef(histaveid, 'phi', 'geopotentiel moyenne', '-', iip1, jjp1, &
61 horiid, llm, 1, llm, zvertiid, 'ave(X)', t_ops, t_wrt)
62
63 ! Traceurs
64 DO iq = 1, nq
65 call histdef(histaveid, ttext(iq), ttext(iq), '-', iip1, jjp1, &
66 horiid, llm, 1, llm, zvertiid, 'ave(X)', t_ops, t_wrt)
67 enddo
68
69 call histdef(histaveid, 'masse', 'masse', 'kg', iip1, jjp1, horiid, 1, &
70 1, 1, -99, 'ave(X)', t_ops, t_wrt)
71 call histdef(histaveid, 'ps', 'pression naturelle au sol', 'Pa', iip1, &
72 jjp1, horiid, 1, 1, 1, -99, 'ave(X)', t_ops, t_wrt)
73 call histdef(histaveid, 'phis', 'geopotentiel au sol', '-', iip1, jjp1, &
74 horiid, 1, 1, 1, -99, 'ave(X)', t_ops, t_wrt)
75
76 call histend(histaveid)
77
78 end subroutine initdynav
79
80 end module initdynav_m

  ViewVC Help
Powered by ViewVC 1.1.21