source: LMDZ5/trunk/libf/phylmd/moy_undefSTD.F90 @ 2313

Last change on this file since 2313 was 2313, checked in by musat, 9 years ago

Suite corrections missing_value pour compilation en mode debug
JG/IM

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 4.9 KB
Line 
1
2! $Id: moy_undefSTD.F90 2313 2015-06-25 12:59:12Z musat $
3
4SUBROUTINE moy_undefstd(itap, itapm1)
5  USE netcdf
6  USE dimphy
7#ifdef CPP_IOIPSL
8  USE phys_state_var_mod
9#endif
10#ifdef CPP_XIOS
11  USE wxios, ONLY: missing_val
12#endif
13 
14  USE phys_cal_mod, ONLY: mth_len
15  IMPLICIT NONE
16  include "clesphys.h"
17#ifndef CPP_XIOS
18  REAL :: missing_val
19#endif
20
21  ! ====================================================================
22
23  ! I. Musat : 09.2004
24
25  ! Moyenne - a des frequences differentes - des valeurs bien definies
26  ! (.NE.missing_val) des variables interpolees a un niveau de
27  ! pression.
28  ! 1) les variables de type "day" (nout=1) ou "mth" (nout=2) sont sommees
29  ! tous les pas de temps de la physique
30
31  ! 2) les variables de type "NMC" (nout=3) sont calculees a partir
32  ! des valeurs instantannees toutes les 6 heures
33
34
35  ! NB: mettre "inst(X)" dans le write_hist*NMC.h !
36  ! ====================================================================
37  ! ym#include "dimensions.h"
38  ! ym      integer jjmp1
39  ! ym      parameter (jjmp1=jjm+1-1/jjm)
40  ! ym#include "dimphy.h"
41
42
43  ! variables Input
44  ! INTEGER nlevSTD, klevSTD, itap
45  ! PARAMETER(klevSTD=17)
46  INTEGER itap, itapm1
47
48  ! variables locales
49  ! INTEGER i, k, nout, n
50  ! PARAMETER(nout=3) !nout=1 day/nout=2 mth/nout=3 NMC
51  INTEGER i, k, n
52  ! REAL dtime, freq_outNMC(nout), freq_moyNMC(nout)
53  ! REAL freq_outNMC(nout), freq_calNMC(nout)
54  REAL freq_moynmc(nout)
55
56  ! variables Output
57  ! REAL tnondef(klon,klevSTD,nout)
58  ! REAL tsumSTD(klon,klevSTD,nout)
59
60  REAL un_jour
61  PARAMETER (un_jour=86400.)
62! REAL missing_val
63
64! missing_val = nf90_fill_real
65#ifndef CPP_XIOS
66      missing_val=missing_val_nf90
67#endif
68
69  DO n = 1, nout
70    IF (freq_outnmc(n)<0) THEN
71      freq_moynmc(n) = (mth_len*un_jour)/freq_calnmc(n)
72      ! print*,'moy_undefSTD n freq_out freq_moy =',
73      ! $n,freq_moyNMC(n)
74    ELSE
75      freq_moynmc(n) = freq_outnmc(n)/freq_calnmc(n)
76    END IF
77
78    ! calcul 1 fois pas mois, 1 fois par jour ou toutes les 6h
79
80    IF (n==1 .AND. itap==itapm1 .OR. n>1 .AND. mod(itap,nint(freq_outnmc(n)/ &
81        dtime))==0) THEN
82
83      ! print*,'moy_undefSTD n itap itapm1',n,itap,itapm1
84
85      DO k = 1, nlevstd
86        DO i = 1, klon
87          IF (tnondef(i,k,n)/=(freq_moynmc(n))) THEN
88            tsumstd(i, k, n) = tsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k,n &
89              ))
90            usumstd(i, k, n) = usumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k,n &
91              ))
92            vsumstd(i, k, n) = vsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k,n &
93              ))
94            wsumstd(i, k, n) = wsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k,n &
95              ))
96            phisumstd(i, k, n) = phisumstd(i, k, n)/ &
97              (freq_moynmc(n)-tnondef(i,k,n))
98            qsumstd(i, k, n) = qsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k,n &
99              ))
100            rhsumstd(i, k, n) = rhsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
101              ,n))
102            uvsumstd(i, k, n) = uvsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
103              ,n))
104            vqsumstd(i, k, n) = vqsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
105              ,n))
106            vtsumstd(i, k, n) = vtsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
107              ,n))
108            wqsumstd(i, k, n) = wqsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
109              ,n))
110            vphisumstd(i, k, n) = vphisumstd(i, k, n)/ &
111              (freq_moynmc(n)-tnondef(i,k,n))
112            wtsumstd(i, k, n) = wtsumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
113              ,n))
114            u2sumstd(i, k, n) = u2sumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
115              ,n))
116            v2sumstd(i, k, n) = v2sumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
117              ,n))
118            t2sumstd(i, k, n) = t2sumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
119              ,n))
120            o3sumstd(i, k, n) = o3sumstd(i, k, n)/(freq_moynmc(n)-tnondef(i,k &
121              ,n))
122            o3daysumstd(i, k, n) = o3daysumstd(i, k, n)/ &
123              (freq_moynmc(n)-tnondef(i,k,n))
124          ELSE
125            tsumstd(i, k, n) = missing_val
126            usumstd(i, k, n) = missing_val
127            vsumstd(i, k, n) = missing_val
128            wsumstd(i, k, n) = missing_val
129            phisumstd(i, k, n) = missing_val
130            qsumstd(i, k, n) = missing_val
131            rhsumstd(i, k, n) = missing_val
132            uvsumstd(i, k, n) = missing_val
133            vqsumstd(i, k, n) = missing_val
134            vtsumstd(i, k, n) = missing_val
135            wqsumstd(i, k, n) = missing_val
136            vphisumstd(i, k, n) = missing_val
137            wtsumstd(i, k, n) = missing_val
138            u2sumstd(i, k, n) = missing_val
139            v2sumstd(i, k, n) = missing_val
140            t2sumstd(i, k, n) = missing_val
141            o3sumstd(i, k, n) = missing_val
142            o3daysumstd(i, k, n) = missing_val
143          END IF !tnondef(i,k,n).NE.(freq_moyNMC(n))
144        END DO !i
145      END DO !k
146    END IF !MOD(itap,NINT(freq_outNMC(n)/dtime)).EQ.0
147
148  END DO !n
149
150  RETURN
151END SUBROUTINE moy_undefstd
Note: See TracBrowser for help on using the repository browser.