source: LMDZ4/branches/V3_test/libf/phylmd/tetalevel.F @ 735

Last change on this file since 735 was 704, checked in by Laurent Fairhead, 18 years ago

Inclusion des modifs de Y. Meurdesoif pour la version V3
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 3.9 KB
Line 
1c================================================================
2c================================================================
3      SUBROUTINE tetalevel(ilon,ilev,lnew,pgcm,pres,Qgcm,Qpres)
4c================================================================
5c================================================================
6      USE dimphy
7      IMPLICIT none
8
9cym#include "dimensions.h"
10cym#include "dimphy.h"
11
12c================================================================
13c
14c Interpoler des champs 3-D u, v et g du modele a un niveau de
15c pression donnee (pres)
16c
17c INPUT:  ilon ----- nombre de points
18c         ilev ----- nombre de couches
19c         lnew ----- true si on doit reinitialiser les poids
20c         pgcm ----- pressions modeles
21c         pres ----- pression vers laquelle on interpolle
22c         Qgcm ----- champ GCM
23c         Qpres ---- champ interpolle au niveau pres
24c
25c================================================================
26c
27c   arguments :
28c   -----------
29
30      INTEGER ilon, ilev
31      logical lnew
32
33      REAL pgcm(ilon,ilev)
34      REAL Qgcm(ilon,ilev)
35      real pres
36      REAL Qpres(ilon)
37
38c   local :
39c   -------
40c
41cym#include "paramet.h"
42c
43      INTEGER,ALLOCATABLE,SAVE :: lt(:), lb(:)
44      REAL,ALLOCATABLE,SAVE    :: aist(:), aisb(:)
45      REAL,SAVE :: ptop, pbot
46      LOGICAL,SAVE :: first = .TRUE.
47c$OMP THREADPRIVATE(lt,lb,aist,aisb,ptop, pbot,first)
48
49      INTEGER i, k
50c
51c     PRINT*,'tetalevel pres=',pres
52      IF (first) THEN
53        ALLOCATE(lt(ilon), lb(ilon))
54        ALLOCATE(aist(ilon), aisb(ilon))
55       
56        first=.FALSE.
57      ENDIF
58c=====================================================================
59      if (lnew) then
60c   on r�nitialise les r�ndicages et les poids
61c=====================================================================
62
63
64c Chercher les 2 couches les plus proches du niveau a obtenir
65c
66c Eventuellement, faire l'extrapolation a partir des deux couches
67c les plus basses ou les deux couches les plus hautes:
68      DO 130 i = 1, ilon
69cIM      IF ( ABS(pres-pgcm(i,ilev) ) .LT.
70         IF ( ABS(pres-pgcm(i,ilev) ) .GT.
71     .        ABS(pres-pgcm(i,1)) ) THEN
72            lt(i) = ilev     ! 2
73            lb(i) = ilev-1   ! 1
74         ELSE
75            lt(i) = 2
76            lb(i) = 1
77         ENDIF
78cIM   PRINT*,'i, ABS(pres-pgcm),ABS(pres-pgcm)',
79cIM  .i, ABS(pres-pgcm(i,ilev)),ABS(pres-pgcm(i,1))
80  130 CONTINUE
81      DO 150 k = 1, ilev-1
82         DO 140 i = 1, ilon
83            pbot = pgcm(i,k)
84            ptop = pgcm(i,k+1)
85cIM         IF (ptop.LE.pres .AND. pbot.GE.pres) THEN
86            IF (ptop.GE.pres .AND. pbot.LE.pres) THEN
87               lt(i) = k+1
88               lb(i) = k
89            ENDIF
90  140    CONTINUE
91  150 CONTINUE
92c
93c Interpolation lineaire:
94c
95      DO i = 1, ilon
96c interpolation en logarithme de pression:
97c
98c ...   Modif . P. Le Van    ( 20/01/98) ....
99c       Modif Fr��ic Hourdin (3/01/02)
100
101c       IF(pgcm(i,lb(i)).NE.0.OR.
102c    $     pgcm(i,lt(i)).NE.0.) THEN
103c
104c       PRINT*,'i,lb,lt,2pgcm,pres',i,lb(i),
105c    .  lt(i),pgcm(i,lb(i)),pgcm(i,lt(i)),pres
106c
107        aist(i) = LOG( pgcm(i,lb(i))/ pres )
108     .       / LOG( pgcm(i,lb(i))/ pgcm(i,lt(i)) )
109        aisb(i) = LOG( pres / pgcm(i,lt(i)) )
110     .       / LOG( pgcm(i,lb(i))/ pgcm(i,lt(i)))
111      enddo
112
113
114      endif ! lnew
115
116c======================================================================
117c    inteprollation
118c======================================================================
119
120      do i=1,ilon
121         Qpres(i)= Qgcm(i,lb(i))*aisb(i)+Qgcm(i,lt(i))*aist(i)
122cIM      PRINT*,'i,Qgcm,Qpres',i,Qgcm(i,lb(i)),aisb(i),
123cIM  $   Qgcm(i,lt(i)),aist(i),Qpres(i)
124      enddo
125c
126c Je mets les vents a zero quand je rencontre une montagne
127      do i = 1, ilon
128cIM      if (pgcm(i,1).LT.pres) THEN
129         if (pgcm(i,1).GT.pres) THEN
130c           Qpres(i)=1e33
131            Qpres(i)=1e+20
132cIM         PRINT*,'i,pgcm(i,1),pres =',i,pgcm(i,1),pres
133         endif
134      enddo
135
136c
137      RETURN
138      END
Note: See TracBrowser for help on using the repository browser.