/[lmdze]/trunk/dyn3d/guide.f
ViewVC logotype

Diff of /trunk/dyn3d/guide.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

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

Legend:
Removed from v.82  
changed lines
  Added in v.91

  ViewVC Help
Powered by ViewVC 1.1.21