source: lmdz_wrf/WRFV3/dyn_nmm/module_PRECIP_ADJUST.F @ 1

Last change on this file since 1 was 1, checked in by lfita, 10 years ago
  • -- --- Opening of the WRF+LMDZ coupling repository --- -- -

WRF: version v3.3
LMDZ: version v1818

More details in:

File size: 10.3 KB
Line 
1!-----------------------------------------------------------------------
2!
3!NCEP_MESO:MODEL_LAYER: PHYSICS
4!
5!----------------------------------------------------------------------
6#include "nmm_loop_basemacros.h"
7#include "nmm_loop_macros.h"
8!-----------------------------------------------------------------------
9!
10      MODULE MODULE_PRECIP_ADJUST
11!
12! This module contains 3 subroutines:
13!     READPCP
14!     CHKSNOW
15!     ADJPPT
16!-----------------------------------------------------------------------
17!***
18!***  Specify the diagnostic point here: (i,j) and the processor number.
19!***  Remember that in WRF, local and global (i,j) are the same, so don't
20!***  use the "local(i,j)" output from glb2loc.f; use the GLOBAL (I,J)
21!***  and the PE_WRF.
22!***
23!
24      INTEGER :: ITEST=346,JTEST=256,TESTPE=53
25!-----------------------------------------------------------------------
26!
27      CONTAINS
28!
29!-----------------------------------------------------------------------
30      SUBROUTINE READPCP(PPTDAT,DDATA,LSPA                              &
31     &  ,IDS,IDE,JDS,JDE,KDS,KDE                                        &
32     &  ,IMS,IME,JMS,JME,KMS,KME                                        &
33     &  ,ITS,ITE,JTS,JTE,KTS,KTE)
34!
35!     ****************************************************************
36!     *                                                              *
37!     *   PRECIPITATION ASSIMILATION INITIALIZATION.                 *
38!     *    READ IN PRECIP ANALYSIS AND DATA MASK AND SET UP ALL      *
39!     *    APPROPRIATE VARIABLES.                                    *
40!     *                   MIKE BALDWIN, MARCH 1994                   *
41!     *                   Adapted to 2-D code, Ying Lin, Mar 1996    *
42!     *                   For WRF/NMM: Y.Lin Mar 2005                *
43!     *                                                              *
44!     ****************************************************************
45!-----------------------------------------------------------------------
46!
47! READ THE BINARY VERSION OF THE PRECIP ANALYSIS.
48!
49      IMPLICIT NONE
50      INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE,                    &
51     &                      IMS,IME,JMS,JME,KMS,KME,                    &
52     &                      ITS,ITE,JTS,JTE,KTS,KTE
53      REAL,DIMENSION(IDS:IDE,JDS:JDE) :: TEMPG
54      REAL,DIMENSION(IMS:IME,JMS:JME) :: TEMPL
55      REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(OUT) :: DDATA, LSPA
56      REAL,DIMENSION(IMS:IME,JMS:JME,3),INTENT(OUT) :: PPTDAT
57      INTEGER :: I, J, IHR
58      INTEGER :: MYPE
59      CHARACTER*256 :: MESSAGE
60!
61! Get the value of MYPE:
62!
63      CALL WRF_GET_MYPROC(MYPE)
64!
65      TEMPG=999.
66!
67      DO IHR=1,3
68        IF(MYPE==0)THEN
69          READ(40+IHR) ((TEMPG(I,J),I=IDS,IDE-1),J=JDS,JDE-1)
70          WRITE(MESSAGE,*) 'IHR=', IHR, ' FINISHED READING PCP TO TEMPG'
71          CALL WRF_MESSAGE(MESSAGE)
72          CLOSE(40+IHR)
73!
74          DO J=JDS,JDE-1
75            DO I=IDS,IDE-1
76! In the binary version of the precip data, missing data are denoted as '999.'
77! Convert the valid data from mm to m:
78              IF (TEMPG(I,J).LT.900.) TEMPG(I,J)=TEMPG(I,J)*0.001
79            ENDDO
80          ENDDO
81        ENDIF
82!
83! Distribute to local temp array:
84        CALL DSTRB(TEMPG,TEMPL,1,1,1,1,1                                &
85     &,                IDS,IDE,JDS,JDE,KDS,KDE                          &
86     &,                IMS,IME,JMS,JME,KMS,KME                          &
87     &,                ITS,ITE,JTS,JTE,KTS,KTE)
88!
89! Place into correct hour slot in PPTDAT:
90        DO J=JMS,JME
91          DO I=IMS,IME
92            PPTDAT(I,J,IHR)=TEMPL(I,J)
93          ENDDO
94        ENDDO
95!
96        IF(MYPE==TESTPE)THEN
97          WRITE(MESSAGE,*) 'ADJPPT-READPCP, IHR',IHR, 'PPTDAT=',        &
98     &      PPTDAT(ITEST,JTEST,IHR)
99          CALL WRF_MESSAGE(MESSAGE)
100        ENDIF
101
102      ENDDO
103!
104! Give DDATA (hourly precipitation analysis partitioned into each physics
105! timestep; partitioning done in ADJPPT) an initial value of 999, because
106! TURBL/SURFCE is called before ADJPPT.  Also initialize LSPA to zero.
107!
108      DDATA=999.
109      LSPA=0.
110!
111      RETURN
112      END SUBROUTINE READPCP
113!
114      SUBROUTINE CHKSNOW(NTSD,DT,NPHS,SR,PPTDAT                         &
115     &  ,IDS,IDE,JDS,JDE,KDS,KDE                                        &
116     &  ,IMS,IME,JMS,JME,KMS,KME                                        &
117     &  ,ITS,ITE,JTS,JTE,KTS,KTE)
118!
119! AT THE FIRST PHYSICS TIME STEP AFTER THE TOP OF EACH HOUR, CHECK THE SNOW
120! ARRAY AGAINST THE SR (SNOW/TOTAL PRECIP RATIO).  IF SR .GE. 0.9, SET THIS
121! POINT TO MISSING (SO WE WON'T DO SNOW ADJUSTMENT HERE).
122!
123!-----------------------------------------------------------------------
124!
125      IMPLICIT NONE
126!
127!-----------------------------------------------------------------------
128!
129      INTEGER,INTENT(IN) :: NTSD,NPHS
130      INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE,                    &
131     &                      IMS,IME,JMS,JME,KMS,KME,                    &
132     &                      ITS,ITE,JTS,JTE,KTS,KTE
133      REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: SR
134      REAL,DIMENSION(IMS:IME,JMS:JME,3),INTENT(INOUT) :: PPTDAT
135      REAL,INTENT(IN) :: DT
136      REAL :: TIMES
137      INTEGER :: I, J, IHR
138      INTEGER :: MYPE
139      CHARACTER*256 :: MESSAGE
140!-----------------------------------------------------------------------
141      TIMES=NTSD*DT
142      IF (MOD(TIMES,3600.) < NPHS*DT) THEN
143        IHR=INT(TIMES)/3600+1
144        IF (IHR > 3) go to 10
145        DO J=MYJS2,MYJE2
146        DO I=MYIS1,MYIE1
147          IF (SR(I,J) >= 0.9) PPTDAT(I,J,IHR) = 999.
148        ENDDO
149        ENDDO
150!
151! Get the value of MYPE:
152!
153        CALL WRF_GET_MYPROC(MYPE)
154!
155        IF (MYPE==TESTPE) THEN
156          WRITE(MESSAGE,1010) TIMES,SR(ITEST,JTEST)
157 1010     FORMAT('ADJPPT-CHKSNOW: TIMES, SR=',F6.0,1X,F6.4)
158          CALL WRF_MESSAGE(MESSAGE)
159        ENDIF
160      ENDIF
161 10   CONTINUE
162      RETURN
163      END SUBROUTINE CHKSNOW
164!
165      SUBROUTINE ADJPPT(NTSD,DT,NPHS,PREC,LSPA,PPTDAT,DDATA             &
166     &  ,IDS,IDE,JDS,JDE,KDS,KDE                                        &
167     &  ,IMS,IME,JMS,JME,KMS,KME                                        &
168     &  ,ITS,ITE,JTS,JTE,KTS,KTE)
169
170!***********************************************************************
171!$$$  SUBPROGRAM DOCUMENTATION BLOCK
172!                .      .    .     
173! SUBPROGRAM:    ADJPPT    PRECIPITATION/CLOUD ADJUSTMENT
174!    PRGRMMR:    Y. LIN    ORG: W/NP22     DATE: 2005/03/30
175!     
176! ABSTRACT:
177!     ADJPPT  MAKES ADJUSTMENT TO MODEL'S TEMPERATURE, MOISTURE, HYDROMETEOR
178!     FIELDS TO BE MORE CONSISTENT WITH THE OBSERVED PRECIPITATION AND CLOUD
179!     TOP PRESSURE
180!     
181!     FOR NOW, AS A FIRST STEP, JUST PARTITION THE INPUT HOURLY PRECIPITATION
182!     OBSERVATION INTO TIME STEPS, AND FEED IT INTO THE SOIL.
183! PROGRAM HISTORY LOG:
184!
185!   2005/03/30  LIN      - BAREBONES PRECIPITATION PARTITION/FEEDING TO GROUND
186! ATTRIBUTES:
187!   LANGUAGE: FORTRAN 90
188!   MACHINE : IBM
189!$$$ 
190!-----------------------------------------------------------------------
191!
192      IMPLICIT NONE
193!
194!-----------------------------------------------------------------------
195      INTEGER,INTENT(IN) :: NPHS, NTSD
196      INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE,                    &
197     &                      IMS,IME,JMS,JME,KMS,KME,                    &
198     &                      ITS,ITE,JTS,JTE,KTS,KTE
199      REAL,INTENT(IN) :: DT
200      REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: PREC
201      REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(OUT) :: DDATA, LSPA
202      REAL,DIMENSION(IMS:IME,JMS:JME,3),INTENT(OUT) :: PPTDAT
203!-----------------------------------------------------------------------
204!***
205!***  LOCAL VARIABLES
206!***
207!-----------------------------------------------------------------------
208      REAL :: DTPHS, FRACT, FRACT1, FRACT2, TIMES, TPHS1, TPHS2
209      INTEGER :: I, J, IHR, IHR1, IHR2, NTSP
210      INTEGER :: MYPE
211      CHARACTER*256 :: MESSAGE
212!
213! Get the value of MYPE:
214!
215      CALL WRF_GET_MYPROC(MYPE)
216!
217      TIMES=NTSD*DT
218      IHR=INT(TIMES)/3600+1
219! Size of physics time step:
220      DTPHS=NPHS*DT
221!
222! Compute the beginning and ending time of the current physics time step,
223! TPHS1 and TPHS2:
224! 
225      NTSP=NTSD/NPHS+1
226      TPHS1=(NTSP-1)*DTPHS
227      TPHS2=NTSP*DTPHS
228!
229      IHR1=INT(TPHS1)/3600+1
230      IHR2=INT(TPHS2)/3600+1
231!
232! Fraction of an hour that falls into IHR1 and IHR2.  Note that IHR1 and IHR2
233! might be identical.
234      IF (IHR1 > 3) THEN
235        GO TO 200
236      ELSEIF (IHR2 > 3) THEN
237        IHR2=3
238        FRACT1=(3600.- MOD(INT(TPHS1),3600))/3600.
239        FRACT2=0.
240      ELSEIF (IHR1 .EQ. IHR2) THEN
241         FRACT1=0.5*DTPHS/3600.
242         FRACT2=FRACT1
243      ELSE
244         FRACT1=(3600.- MOD(INT(TPHS1),3600))/3600.
245         FRACT2=FLOAT(MOD(INT(TPHS2),3600))/3600.
246      ENDIF
247!
248      FRACT=FRACT1 + FRACT2
249!
250      IF (MYPE==TESTPE) THEN
251         WRITE(MESSAGE,1010) NTSD,NTSP,TIMES,IHR1,IHR2,TPHS1,TPHS2,      &
252      &    FRACT1,FRACT2
253 1010    FORMAT('ADJPPT: NTSD,NTSP,TIMES=',I4,1X,I4,1X,F6.0,' IHR1,IHR2=' &
254      &   ,I1,1X,I1,' TPHS1,TPHS2=',F6.0,1X,F6.0,' FRACT1,FRACT2='        &
255      &   ,2(1X,F6.4))
256        CALL WRF_MESSAGE(MESSAGE)
257      ENDIF
258!
259!-----------------------------------------------------------------------
260!   FRACT1/2 IS THE FRACTION OF IHR1/2'S PRECIP THAT WE WANT FOR
261!   THIS ADJUSTMENT (assuming that the physics time step spans over
262!   IHR1 and IHR2.  If not, then IHR1=IHR2).
263!-----------------------------------------------------------------------
264!   SET UP OBSERVED PRECIP FOR THIS TIMESTEP IN DDATA
265!-----------------------------------------------------------------------
266      DO J=MYJS2,MYJE2
267      DO I=MYIS1,MYIE1
268! Note sometimes IHR1=IHR2. 
269        IF (PPTDAT(I,J,IHR1).GT.900..OR.PPTDAT(I,J,IHR2).GT.900.) THEN
270          DDATA(I,J) = 999.
271          LSPA(I,J) = LSPA(I,J) + PREC(I,J)
272          GO TO 100
273        ELSE
274          IF (IHR2 .LE. 3) then
275            DDATA(I,J) = PPTDAT(I,J,IHR1)*FRACT1                        &
276     &                 + PPTDAT(I,J,IHR2)*FRACT2
277          ELSE
278            DDATA(I,J) = PPTDAT(I,J,IHR1)*FRACT1
279          ENDIF
280!
281           LSPA(I,J) = LSPA(I,J) + DDATA(I,J)
282        ENDIF
283        IF (I.EQ.ITEST .AND. J.EQ.JTEST .AND. MYPE.EQ.TESTPE) THEN
284          WRITE(MESSAGE,1020) DDATA(I,J), PREC(I,J), LSPA(I,J)
285 1020     FORMAT('ADJPPT: DDATA=',E12.6, ' PREC=',E12.6,' LSPA=',E12.6)
286          CALL WRF_MESSAGE(MESSAGE)
287        ENDIF
288!
289 100    CONTINUE
290      ENDDO
291      ENDDO
292!
293 200  CONTINUE
294
295      RETURN
296      END SUBROUTINE ADJPPT
297END MODULE module_PRECIP_ADJUST
Note: See TracBrowser for help on using the repository browser.