trunk/dyn3d/guide.f

MODULE guide_m

  ! From dyn3d/guide.F, version 1.3 2005/05/25 13:10:09
  ! and dyn3d/guide.h, version 1.1.1.1 2004/05/19 12:53:06

  IMPLICIT NONE

  REAL aire_min, aire_max

CONTAINS

  SUBROUTINE guide(itau, ucov, vcov, teta, q, masse, ps)

    ! Author: F.Hourdin

    USE comconst, ONLY : cpp, daysec, dtvr, kappa
    USE comgeom, ONLY : aire, rlatu, rlonv
    USE comvert, ONLY : ap, bp, preff, presnivs
    USE conf_gcm_m, ONLY : day_step, iperiod
    use conf_guide_m, only: conf_guide, guide_u, guide_v, guide_t, guide_q, &
         guide_p, ncep, ini_anal, tau_min_u, tau_max_u, tau_min_v, tau_max_v, &
         tau_min_t, tau_max_t, tau_min_q, tau_max_q, tau_min_p, tau_max_p, &
         online
    USE dimens_m, ONLY : jjm, llm
    USE exner_hyb_m, ONLY : exner_hyb
    USE inigrads_m, ONLY : inigrads
    use netcdf, only: nf90_nowrite, nf90_open, nf90_close, nf90_inq_dimid, &
         nf90_inquire_dimension
    use nr_util, only: pi
    USE paramet_m, ONLY : iip1, ip1jm, ip1jmp1, jjp1, llmp1
    USE q_sat_m, ONLY : q_sat
    USE serre, ONLY : clat, clon
    use tau2alpha_m, only: tau2alpha, dxdys

    ! variables dynamiques
    REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm) ! vents covariants
    REAL, intent(inout):: teta(ip1jmp1, llm) ! temperature potentielle 
    REAL q(ip1jmp1, llm) ! temperature potentielle 
    REAL ps(ip1jmp1) ! pression au sol
    REAL masse(ip1jmp1, llm) ! masse d'air

    ! variables dynamiques pour les reanalyses.
    REAL, save:: ucovrea1(ip1jmp1, llm), vcovrea1(ip1jm, llm) !vts cov reas
    REAL, save:: tetarea1(ip1jmp1, llm) ! temp pot reales
    REAL, save:: qrea1(ip1jmp1, llm) ! temp pot reales
    REAL, save:: psrea1(ip1jmp1) ! ps
    REAL, save:: ucovrea2(ip1jmp1, llm), vcovrea2(ip1jm, llm) !vts cov reas
    REAL, save:: tetarea2(ip1jmp1, llm) ! temp pot reales
    REAL, save:: qrea2(ip1jmp1, llm) ! temp pot reales
    REAL, save:: masserea2(ip1jmp1, llm) ! masse
    REAL, save:: psrea2(ip1jmp1) ! ps

    REAL, save:: alpha_q(ip1jmp1)
    REAL, save:: alpha_t(ip1jmp1), alpha_p(ip1jmp1)
    REAL, save:: alpha_u(ip1jmp1), alpha_v(ip1jm)
    REAL dday_step, toto, reste
    real, save:: itau_test
    INTEGER, save:: step_rea, count_no_rea

    INTEGER ilon, ilat
    REAL factt, ztau(ip1jmp1)

    INTEGER, INTENT(IN):: itau
    INTEGER ij, l
    INTEGER ncidpl, varidpl, status
    INTEGER rcod, rid
    REAL ditau, tau, a
    INTEGER, SAVE:: nlev

    ! TEST SUR QSAT
    REAL p(ip1jmp1, llmp1), pk(ip1jmp1, llm), pks(ip1jmp1)
    REAL pkf(ip1jmp1, llm)
    REAL pres(ip1jmp1, llm)

    REAL qsat(ip1jmp1, llm)
    REAL unskap
    REAL tnat(ip1jmp1, llm)

    LOGICAL:: first = .TRUE.
    CHARACTER(len=10) file
    INTEGER:: igrads = 2
    REAL:: dtgrads = 100.

    !-----------------------------------------------------------------------

    PRINT *, 'Call sequence information: guide'

    ! calcul de l'humidite saturante

    forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
    CALL massdair(p, masse)
    CALL exner_hyb(ps, p, pks, pk, pkf)
    tnat(:, :) = pk(:, :)*teta(:, :)/cpp
    unskap = 1./kappa
    pres(:, :) = preff*(pk(:, :)/cpp)**unskap
    qsat = q_sat(tnat, pres)

    ! initialisations pour la lecture des reanalyses.
    ! alpha determine la part des injections de donnees a chaque etape
    ! alpha=1 signifie pas d'injection
    ! alpha=0 signifie injection totale

    IF (online==-1) THEN
       RETURN
    END IF

    IF (first) THEN
       CALL conf_guide
       file = 'guide'
       CALL inigrads(igrads, rlonv, 180./pi, -180., 180., rlatu, -90., 90., &
            180./pi, presnivs, 1., dtgrads, file, 'dyn_zon ')
       PRINT *, '1: en-ligne, 0: hors-ligne (x=x_rea), -1: climat (x=x_gcm)'
       IF (online==-1) RETURN

       IF (online==1) THEN
          ! Constantes de temps de rappel en jour
          ! 0.1 c'est en gros 2h30. 
          ! 1e10 est une constante infinie donc en gros pas de guidage

          ! coordonnees du centre du zoom
          CALL coordij(clon, clat, ilon, ilat)
          ! aire de la maille au centre du zoom
          aire_min = aire(ilon+(ilat-1)*iip1)
          ! aire maximale de la maille
          aire_max = 0.
          DO ij = 1, ip1jmp1
             aire_max = max(aire_max, aire(ij))
          END DO
          ! factt = pas de temps en fraction de jour
          factt = dtvr*iperiod/daysec

          CALL tau2alpha(3, iip1, jjm, factt, tau_min_v, tau_max_v, alpha_v)
          CALL tau2alpha(2, iip1, jjp1, factt, tau_min_u, tau_max_u, alpha_u)
          CALL tau2alpha(1, iip1, jjp1, factt, tau_min_t, tau_max_t, alpha_t)
          CALL tau2alpha(1, iip1, jjp1, factt, tau_min_p, tau_max_p, alpha_p)
          CALL tau2alpha(1, iip1, jjp1, factt, tau_min_q, tau_max_q, alpha_q)

          CALL dump2d(iip1, jjp1, aire, 'AIRE MAILLe ')
          CALL dump2d(iip1, jjp1, alpha_u, 'COEFF U ')
          CALL dump2d(iip1, jjp1, alpha_t, 'COEFF T ')

          ! Cas ou on force exactement par les variables analysees
       ELSE
          alpha_t = 0.
          alpha_u = 0.
          alpha_v = 0.
          alpha_p = 0.
          ! physic=.false.
       END IF

       itau_test = 1001
       step_rea = 1
       count_no_rea = 0
       ncidpl = -99

       ! itau_test montre si l'importation a deja ete faite au rang itau
       ! lecture d'un fichier netcdf pour determiner le nombre de niveaux
       if (guide_u) then
          if (ncidpl.eq.-99) rcod=nf90_open('u.nc',Nf90_NOWRITe,ncidpl)
       endif

       if (guide_v) then
          if (ncidpl.eq.-99) rcod=nf90_open('v.nc',nf90_nowrite,ncidpl)
       endif

       if (guide_T) then
          if (ncidpl.eq.-99) rcod=nf90_open('T.nc',nf90_nowrite,ncidpl)
       endif

       if (guide_Q) then
          if (ncidpl.eq.-99) rcod=nf90_open('hur.nc',nf90_nowrite, ncidpl)
       endif

       IF (ncep) THEN
          status = nf90_inq_dimid(ncidpl, 'LEVEL', rid)
       ELSE
          status = nf90_inq_dimid(ncidpl, 'PRESSURE', rid)
       END IF
       status = nf90_inquire_dimension(ncidpl, rid, len=nlev)
       PRINT *, 'nlev', nlev
       rcod = nf90_close(ncidpl)
       ! Lecture du premier etat des reanalyses.
       CALL read_reanalyse(1, ps, ucovrea2, vcovrea2, tetarea2, qrea2, &
            masserea2, psrea2, 1, nlev)
       qrea2(:, :) = max(qrea2(:, :), 0.1)

       ! Debut de l'integration temporelle:
    END IF ! first

    ! IMPORTATION DES VENTS, PRESSION ET TEMPERATURE REELS:

    ditau = real(itau)
    dday_step = real(day_step)
    WRITE (*, *) 'ditau, dday_step'
    WRITE (*, *) ditau, dday_step
    toto = 4*ditau/dday_step
    reste = toto - aint(toto)

    IF (reste==0.) THEN
       IF (itau_test==itau) THEN
          WRITE (*, *) 'deuxieme passage de advreel a itau=', itau
          STOP
       ELSE
          vcovrea1(:, :) = vcovrea2(:, :)
          ucovrea1(:, :) = ucovrea2(:, :)
          tetarea1(:, :) = tetarea2(:, :)
          qrea1(:, :) = qrea2(:, :)

          PRINT *, 'LECTURE REANALYSES, pas ', step_rea, 'apres ', &
               count_no_rea, ' non lectures'
          step_rea = step_rea + 1
          itau_test = itau
          CALL read_reanalyse(step_rea, ps, ucovrea2, vcovrea2, tetarea2, &
               qrea2, masserea2, psrea2, 1, nlev)
          qrea2(:, :) = max(qrea2(:, :), 0.1)
          factt = dtvr*iperiod/daysec
          ztau(:) = factt/max(alpha_t(:), 1.E-10)
          CALL wrgrads(igrads, 1, aire, 'aire ', 'aire ')
          CALL wrgrads(igrads, 1, dxdys, 'dxdy ', 'dxdy ')
          CALL wrgrads(igrads, 1, alpha_u, 'au ', 'au ')
          CALL wrgrads(igrads, 1, alpha_t, 'at ', 'at ')
          CALL wrgrads(igrads, 1, ztau, 'taut ', 'taut ')
          CALL wrgrads(igrads, llm, ucov, 'u ', 'u ')
          CALL wrgrads(igrads, llm, ucovrea2, 'ua ', 'ua ')
          CALL wrgrads(igrads, llm, teta, 'T ', 'T ')
          CALL wrgrads(igrads, llm, tetarea2, 'Ta ', 'Ta ')
          CALL wrgrads(igrads, llm, qrea2, 'Qa ', 'Qa ')
          CALL wrgrads(igrads, llm, q, 'Q ', 'Q ')

          CALL wrgrads(igrads, llm, qsat, 'QSAT ', 'QSAT ')

       END IF
    ELSE
       count_no_rea = count_no_rea + 1
    END IF

    ! Guidage
    ! x_gcm = a * x_gcm + (1-a) * x_reanalyses

    IF (ini_anal) PRINT *, 'ATTENTION !!! ON PART DU GUIDAGE'

    ditau = real(itau)
    dday_step = real(day_step)

    tau = 4*ditau/dday_step
    tau = tau - aint(tau)

    ! ucov
    IF (guide_u) THEN
       DO l = 1, llm
          DO ij = 1, ip1jmp1
             a = (1.-tau)*ucovrea1(ij, l) + tau*ucovrea2(ij, l)
             ucov(ij, l) = (1.-alpha_u(ij))*ucov(ij, l) + alpha_u(ij)*a
             IF (first .AND. ini_anal) ucov(ij, l) = a
          END DO
       END DO
    END IF

    IF (guide_t) THEN
       DO l = 1, llm
          DO ij = 1, ip1jmp1
             a = (1.-tau)*tetarea1(ij, l) + tau*tetarea2(ij, l)
             teta(ij, l) = (1.-alpha_t(ij))*teta(ij, l) + alpha_t(ij)*a
             IF (first .AND. ini_anal) teta(ij, l) = a
          END DO
       END DO
    END IF

    ! P
    IF (guide_p) THEN
       DO ij = 1, ip1jmp1
          a = (1.-tau)*psrea1(ij) + tau*psrea2(ij)
          ps(ij) = (1.-alpha_p(ij))*ps(ij) + alpha_p(ij)*a
          IF (first .AND. ini_anal) ps(ij) = a
       END DO
       forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
       CALL massdair(p, masse)
    END IF

    ! q
    IF (guide_q) THEN
       DO l = 1, llm
          DO ij = 1, ip1jmp1
             a = (1.-tau)*qrea1(ij, l) + tau*qrea2(ij, l)
             ! hum relative en % -> hum specif
             a = qsat(ij, l)*a*0.01
             q(ij, l) = (1.-alpha_q(ij))*q(ij, l) + alpha_q(ij)*a
             IF (first .AND. ini_anal) q(ij, l) = a
          END DO
       END DO
    END IF

    ! vcov
    IF (guide_v) THEN
       DO l = 1, llm
          DO ij = 1, ip1jm
             a = (1.-tau)*vcovrea1(ij, l) + tau*vcovrea2(ij, l)
             vcov(ij, l) = (1.-alpha_v(ij))*vcov(ij, l) + alpha_v(ij)*a
             IF (first .AND. ini_anal) vcov(ij, l) = a
          END DO
          IF (first .AND. ini_anal) vcov(ij, l) = a
       END DO
    END IF

    first = .FALSE.

  END SUBROUTINE guide

END MODULE guide_m
1	guez	20	MODULE guide_m
2	guez	3
3	guez	29	! From dyn3d/guide.F, version 1.3 2005/05/25 13:10:09
4			! and dyn3d/guide.h, version 1.1.1.1 2004/05/19 12:53:06
5	guez	3
6	guez	37	IMPLICIT NONE
7
8	guez	36	REAL aire_min, aire_max
9	guez	3
10	guez	20	CONTAINS
11	guez	3
12	guez	29	SUBROUTINE guide(itau, ucov, vcov, teta, q, masse, ps)
13	guez	3
14	guez	29	! Author: F.Hourdin
15	guez	3
16	guez	39	USE comconst, ONLY : cpp, daysec, dtvr, kappa
17	guez	44	USE comgeom, ONLY : aire, rlatu, rlonv
18	guez	29	USE comvert, ONLY : ap, bp, preff, presnivs
19			USE conf_gcm_m, ONLY : day_step, iperiod
20	guez	44	use conf_guide_m, only: conf_guide, guide_u, guide_v, guide_t, guide_q, &
21			guide_p, ncep, ini_anal, tau_min_u, tau_max_u, tau_min_v, tau_max_v, &
22			tau_min_t, tau_max_t, tau_min_q, tau_max_q, tau_min_p, tau_max_p, &
23			online
24	guez	39	USE dimens_m, ONLY : jjm, llm
25	guez	29	USE exner_hyb_m, ONLY : exner_hyb
26			USE inigrads_m, ONLY : inigrads
27	guez	44	use netcdf, only: nf90_nowrite, nf90_open, nf90_close, nf90_inq_dimid, &
28			nf90_inquire_dimension
29	guez	39	use nr_util, only: pi
30			USE paramet_m, ONLY : iip1, ip1jm, ip1jmp1, jjp1, llmp1
31			USE q_sat_m, ONLY : q_sat
32			USE serre, ONLY : clat, clon
33	guez	44	use tau2alpha_m, only: tau2alpha, dxdys
34	guez	3
35	guez	44	! variables dynamiques
36	guez	36	REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm) ! vents covariants
37	guez	29	REAL, intent(inout):: teta(ip1jmp1, llm) ! temperature potentielle
38	guez	36	REAL q(ip1jmp1, llm) ! temperature potentielle
39	guez	44	REAL ps(ip1jmp1) ! pression au sol
40	guez	36	REAL masse(ip1jmp1, llm) ! masse d'air
41	guez	3
42	guez	44	! variables dynamiques pour les reanalyses.
43			REAL, save:: ucovrea1(ip1jmp1, llm), vcovrea1(ip1jm, llm) !vts cov reas
44			REAL, save:: tetarea1(ip1jmp1, llm) ! temp pot reales
45			REAL, save:: qrea1(ip1jmp1, llm) ! temp pot reales
46			REAL, save:: psrea1(ip1jmp1) ! ps
47			REAL, save:: ucovrea2(ip1jmp1, llm), vcovrea2(ip1jm, llm) !vts cov reas
48			REAL, save:: tetarea2(ip1jmp1, llm) ! temp pot reales
49			REAL, save:: qrea2(ip1jmp1, llm) ! temp pot reales
50			REAL, save:: masserea2(ip1jmp1, llm) ! masse
51			REAL, save:: psrea2(ip1jmp1) ! ps
52	guez	3
53	guez	44	REAL, save:: alpha_q(ip1jmp1)
54			REAL, save:: alpha_t(ip1jmp1), alpha_p(ip1jmp1)
55			REAL, save:: alpha_u(ip1jmp1), alpha_v(ip1jm)
56			REAL dday_step, toto, reste
57			real, save:: itau_test
58			INTEGER, save:: step_rea, count_no_rea
59	guez	3
60	guez	36	INTEGER ilon, ilat
61			REAL factt, ztau(ip1jmp1)
62	guez	3
63	guez	44	INTEGER, INTENT(IN):: itau
64	guez	36	INTEGER ij, l
65	guez	44	INTEGER ncidpl, varidpl, status
66	guez	36	INTEGER rcod, rid
67			REAL ditau, tau, a
68	guez	44	INTEGER, SAVE:: nlev
69	guez	3
70	guez	44	! TEST SUR QSAT
71	guez	36	REAL p(ip1jmp1, llmp1), pk(ip1jmp1, llm), pks(ip1jmp1)
72			REAL pkf(ip1jmp1, llm)
73			REAL pres(ip1jmp1, llm)
74	guez	3
75	guez	36	REAL qsat(ip1jmp1, llm)
76			REAL unskap
77			REAL tnat(ip1jmp1, llm)
78	guez	3
79	guez	37	LOGICAL:: first = .TRUE.
80	guez	44	CHARACTER(len=10) file
81			INTEGER:: igrads = 2
82			REAL:: dtgrads = 100.
83	guez	3
84	guez	29	!-----------------------------------------------------------------------
85	guez	3
86	guez	29	PRINT *, 'Call sequence information: guide'
87	guez	3
88	guez	29	! calcul de l'humidite saturante
89	guez	3
90	guez	37	forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
91	guez	29	CALL massdair(p, masse)
92			CALL exner_hyb(ps, p, pks, pk, pkf)
93			tnat(:, :) = pk(:, :)*teta(:, :)/cpp
94			unskap = 1./kappa
95			pres(:, :) = preff(pk(:, :)/cpp)*unskap
96			qsat = q_sat(tnat, pres)
97	guez	3
98	guez	44	! initialisations pour la lecture des reanalyses.
99			! alpha determine la part des injections de donnees a chaque etape
100			! alpha=1 signifie pas d'injection
101			! alpha=0 signifie injection totale
102	guez	3
103	guez	29	IF (online==-1) THEN
104			RETURN
105			END IF
106	guez	3
107	guez	29	IF (first) THEN
108			CALL conf_guide
109			file = 'guide'
110			CALL inigrads(igrads, rlonv, 180./pi, -180., 180., rlatu, -90., 90., &
111			180./pi, presnivs, 1., dtgrads, file, 'dyn_zon ')
112	guez	44	PRINT *, '1: en-ligne, 0: hors-ligne (x=x_rea), -1: climat (x=x_gcm)'
113			IF (online==-1) RETURN
114	guez	3
115	guez	29	IF (online==1) THEN
116	guez	44	! Constantes de temps de rappel en jour
117			! 0.1 c'est en gros 2h30.
118			! 1e10 est une constante infinie donc en gros pas de guidage
119	guez	3
120	guez	44	! coordonnees du centre du zoom
121	guez	29	CALL coordij(clon, clat, ilon, ilat)
122	guez	44	! aire de la maille au centre du zoom
123	guez	29	aire_min = aire(ilon+(ilat-1)*iip1)
124	guez	44	! aire maximale de la maille
125	guez	29	aire_max = 0.
126			DO ij = 1, ip1jmp1
127			aire_max = max(aire_max, aire(ij))
128			END DO
129	guez	44	! factt = pas de temps en fraction de jour
130	guez	29	factt = dtvr*iperiod/daysec
131	guez	3
132	guez	29	CALL tau2alpha(3, iip1, jjm, factt, tau_min_v, tau_max_v, alpha_v)
133			CALL tau2alpha(2, iip1, jjp1, factt, tau_min_u, tau_max_u, alpha_u)
134			CALL tau2alpha(1, iip1, jjp1, factt, tau_min_t, tau_max_t, alpha_t)
135			CALL tau2alpha(1, iip1, jjp1, factt, tau_min_p, tau_max_p, alpha_p)
136			CALL tau2alpha(1, iip1, jjp1, factt, tau_min_q, tau_max_q, alpha_q)
137	guez	3
138	guez	29	CALL dump2d(iip1, jjp1, aire, 'AIRE MAILLe ')
139	guez	44	CALL dump2d(iip1, jjp1, alpha_u, 'COEFF U ')
140			CALL dump2d(iip1, jjp1, alpha_t, 'COEFF T ')
141	guez	3
142	guez	44	! Cas ou on force exactement par les variables analysees
143	guez	29	ELSE
144			alpha_t = 0.
145			alpha_u = 0.
146			alpha_v = 0.
147			alpha_p = 0.
148	guez	44	! physic=.false.
149	guez	29	END IF
150	guez	3
151	guez	29	itau_test = 1001
152			step_rea = 1
153			count_no_rea = 0
154			ncidpl = -99
155	guez	3
156	guez	44	! itau_test montre si l'importation a deja ete faite au rang itau
157	guez	29	! lecture d'un fichier netcdf pour determiner le nombre de niveaux
158			if (guide_u) then
159			if (ncidpl.eq.-99) rcod=nf90_open('u.nc',Nf90_NOWRITe,ncidpl)
160			endif
161	guez	3
162	guez	29	if (guide_v) then
163			if (ncidpl.eq.-99) rcod=nf90_open('v.nc',nf90_nowrite,ncidpl)
164			endif
165	guez	3
166	guez	29	if (guide_T) then
167			if (ncidpl.eq.-99) rcod=nf90_open('T.nc',nf90_nowrite,ncidpl)
168			endif
169	guez	3
170	guez	29	if (guide_Q) then
171			if (ncidpl.eq.-99) rcod=nf90_open('hur.nc',nf90_nowrite, ncidpl)
172			endif
173	guez	3
174	guez	29	IF (ncep) THEN
175	guez	44	status = nf90_inq_dimid(ncidpl, 'LEVEL', rid)
176	guez	29	ELSE
177	guez	44	status = nf90_inq_dimid(ncidpl, 'PRESSURE', rid)
178	guez	29	END IF
179	guez	44	status = nf90_inquire_dimension(ncidpl, rid, len=nlev)
180	guez	29	PRINT *, 'nlev', nlev
181			rcod = nf90_close(ncidpl)
182	guez	44	! Lecture du premier etat des reanalyses.
183	guez	29	CALL read_reanalyse(1, ps, ucovrea2, vcovrea2, tetarea2, qrea2, &
184			masserea2, psrea2, 1, nlev)
185			qrea2(:, :) = max(qrea2(:, :), 0.1)
186	guez	3
187	guez	44	! Debut de l'integration temporelle:
188	guez	29	END IF ! first
189	guez	3
190	guez	29	! IMPORTATION DES VENTS, PRESSION ET TEMPERATURE REELS:
191	guez	3
192	guez	29	ditau = real(itau)
193			dday_step = real(day_step)
194			WRITE (, ) 'ditau, dday_step'
195			WRITE (, ) ditau, dday_step
196			toto = 4*ditau/dday_step
197			reste = toto - aint(toto)
198	guez	3
199	guez	29	IF (reste==0.) THEN
200			IF (itau_test==itau) THEN
201			WRITE (, ) 'deuxieme passage de advreel a itau=', itau
202			STOP
203			ELSE
204			vcovrea1(:, :) = vcovrea2(:, :)
205			ucovrea1(:, :) = ucovrea2(:, :)
206			tetarea1(:, :) = tetarea2(:, :)
207			qrea1(:, :) = qrea2(:, :)
208	guez	3
209	guez	29	PRINT *, 'LECTURE REANALYSES, pas ', step_rea, 'apres ', &
210			count_no_rea, ' non lectures'
211			step_rea = step_rea + 1
212			itau_test = itau
213			CALL read_reanalyse(step_rea, ps, ucovrea2, vcovrea2, tetarea2, &
214			qrea2, masserea2, psrea2, 1, nlev)
215			qrea2(:, :) = max(qrea2(:, :), 0.1)
216			factt = dtvr*iperiod/daysec
217			ztau(:) = factt/max(alpha_t(:), 1.E-10)
218	guez	44	CALL wrgrads(igrads, 1, aire, 'aire ', 'aire ')
219			CALL wrgrads(igrads, 1, dxdys, 'dxdy ', 'dxdy ')
220			CALL wrgrads(igrads, 1, alpha_u, 'au ', 'au ')
221			CALL wrgrads(igrads, 1, alpha_t, 'at ', 'at ')
222			CALL wrgrads(igrads, 1, ztau, 'taut ', 'taut ')
223			CALL wrgrads(igrads, llm, ucov, 'u ', 'u ')
224			CALL wrgrads(igrads, llm, ucovrea2, 'ua ', 'ua ')
225			CALL wrgrads(igrads, llm, teta, 'T ', 'T ')
226			CALL wrgrads(igrads, llm, tetarea2, 'Ta ', 'Ta ')
227			CALL wrgrads(igrads, llm, qrea2, 'Qa ', 'Qa ')
228			CALL wrgrads(igrads, llm, q, 'Q ', 'Q ')
229	guez	3
230	guez	44	CALL wrgrads(igrads, llm, qsat, 'QSAT ', 'QSAT ')
231	guez	3
232	guez	29	END IF
233			ELSE
234			count_no_rea = count_no_rea + 1
235			END IF
236	guez	3
237	guez	44	! Guidage
238			! x_gcm = a * x_gcm + (1-a) * x_reanalyses
239	guez	3
240	guez	29	IF (ini_anal) PRINT *, 'ATTENTION !!! ON PART DU GUIDAGE'
241	guez	3
242	guez	29	ditau = real(itau)
243			dday_step = real(day_step)
244	guez	3
245	guez	29	tau = 4*ditau/dday_step
246			tau = tau - aint(tau)
247	guez	3
248	guez	44	! ucov
249	guez	29	IF (guide_u) THEN
250			DO l = 1, llm
251			DO ij = 1, ip1jmp1
252			a = (1.-tau)ucovrea1(ij, l) + tauucovrea2(ij, l)
253			ucov(ij, l) = (1.-alpha_u(ij))ucov(ij, l) + alpha_u(ij)a
254			IF (first .AND. ini_anal) ucov(ij, l) = a
255			END DO
256			END DO
257			END IF
258	guez	3
259	guez	29	IF (guide_t) THEN
260			DO l = 1, llm
261			DO ij = 1, ip1jmp1
262			a = (1.-tau)tetarea1(ij, l) + tautetarea2(ij, l)
263			teta(ij, l) = (1.-alpha_t(ij))teta(ij, l) + alpha_t(ij)a
264			IF (first .AND. ini_anal) teta(ij, l) = a
265			END DO
266			END DO
267			END IF
268	guez	3
269	guez	44	! P
270	guez	29	IF (guide_p) THEN
271			DO ij = 1, ip1jmp1
272			a = (1.-tau)psrea1(ij) + taupsrea2(ij)
273			ps(ij) = (1.-alpha_p(ij))ps(ij) + alpha_p(ij)a
274			IF (first .AND. ini_anal) ps(ij) = a
275			END DO
276	guez	37	forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
277	guez	29	CALL massdair(p, masse)
278			END IF
279	guez	3
280	guez	44	! q
281	guez	29	IF (guide_q) THEN
282			DO l = 1, llm
283			DO ij = 1, ip1jmp1
284			a = (1.-tau)qrea1(ij, l) + tauqrea2(ij, l)
285	guez	44	! hum relative en % -> hum specif
286	guez	29	a = qsat(ij, l)a0.01
287			q(ij, l) = (1.-alpha_q(ij))q(ij, l) + alpha_q(ij)a
288			IF (first .AND. ini_anal) q(ij, l) = a
289			END DO
290			END DO
291			END IF
292	guez	3
293	guez	29	! vcov
294			IF (guide_v) THEN
295			DO l = 1, llm
296			DO ij = 1, ip1jm
297			a = (1.-tau)vcovrea1(ij, l) + tauvcovrea2(ij, l)
298			vcov(ij, l) = (1.-alpha_v(ij))vcov(ij, l) + alpha_v(ij)a
299			IF (first .AND. ini_anal) vcov(ij, l) = a
300			END DO
301			IF (first .AND. ini_anal) vcov(ij, l) = a
302			END DO
303			END IF
304	guez	3
305	guez	29	first = .FALSE.
306	guez	3
307	guez	29	END SUBROUTINE guide
308	guez	20
309			END MODULE guide_m