source: LMDZ5/branches/LMDZ5_SPLA/libf/phylmd/plevel_new.F90 @ 4934

Last change on this file since 4934 was 1992, checked in by lguez, 11 years ago

Converted to free source form files in libf/phylmd which were still in
fixed source form. The conversion was done using the polish mode of
the NAG Fortran Compiler.

In addition to converting to free source form, the processing of the
files also:

-- indented the code (including comments);

-- set Fortran keywords to uppercase, and set all other identifiers
to lower case;

-- added qualifiers to end statements (for example "end subroutine
conflx", instead of "end");

-- changed the terminating statements of all DO loops so that each
loop ends with an ENDDO statement (instead of a labeled continue).

-- replaced #include by include.

  • 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
File size: 3.8 KB
Line 
1
2! $Header: /home/cvsroot/LMDZ4/libf/phylmd/plevel.F,v 1.1.1.1.10.1 2006/08/17
3! 15:41:51 fairhead Exp $
4
5! ================================================================
6! ================================================================
7SUBROUTINE plevel_new(ilon, ilev, klevstd, lnew, pgcm, pres, qgcm, qpres)
8  ! ================================================================
9  ! ================================================================
10  USE netcdf
11  USE dimphy
12  IMPLICIT NONE
13
14  ! ym#include "dimensions.h"
15  ! y#include "dimphy.h"
16
17  ! ================================================================
18
19  ! Interpoler des champs 3-D u, v et g du modele a un niveau de
20  ! pression donnee (pres)
21
22  ! INPUT:  ilon ----- nombre de points
23  ! ilev ----- nombre de couches
24  ! lnew ----- true si on doit reinitialiser les poids
25  ! pgcm ----- pressions modeles
26  ! pres ----- pression vers laquelle on interpolle
27  ! Qgcm ----- champ GCM
28  ! Qpres ---- champ interpolle au niveau pres
29
30  ! ================================================================
31
32  ! arguments :
33  ! -----------
34
35  INTEGER ilon, ilev, klevstd
36  LOGICAL lnew
37
38  REAL pgcm(ilon, ilev)
39  REAL qgcm(ilon, ilev)
40  REAL pres(klevstd)
41  REAL qpres(ilon, klevstd)
42
43  ! local :
44  ! -------
45
46  ! ym      INTEGER lt(klon), lb(klon)
47  ! ym      REAL ptop, pbot, aist(klon), aisb(klon)
48
49  ! ym      save lt,lb,ptop,pbot,aist,aisb
50  INTEGER, ALLOCATABLE, SAVE, DIMENSION (:, :) :: lt, lb
51  REAL, ALLOCATABLE, SAVE, DIMENSION (:, :) :: aist, aisb
52  !$OMP THREADPRIVATE(lt,lb,aist,aisb)
53  REAL, SAVE :: ptop, pbot
54  !$OMP THREADPRIVATE(ptop, pbot)
55  LOGICAL, SAVE :: first = .TRUE.
56  INTEGER :: nlev
57  !$OMP THREADPRIVATE(first)
58  INTEGER i, k
59
60  REAL missing_val
61
62  missing_val = nf90_fill_real
63
64  IF (first) THEN
65    ALLOCATE (lt(klon,klevstd), lb(klon,klevstd))
66    ALLOCATE (aist(klon,klevstd), aisb(klon,klevstd))
67    first = .FALSE.
68  END IF
69
70  ! =====================================================================
71  IF (lnew) THEN
72    ! on reinitialise les reindicages et les poids
73    ! =====================================================================
74
75
76    ! Chercher les 2 couches les plus proches du niveau a obtenir
77
78    ! Eventuellement, faire l'extrapolation a partir des deux couches
79    ! les plus basses ou les deux couches les plus hautes:
80
81
82    DO nlev = 1, klevstd
83      DO i = 1, klon
84        IF (abs(pres(nlev)-pgcm(i,ilev))<abs(pres(nlev)-pgcm(i,1))) THEN
85          lt(i, nlev) = ilev ! 2
86          lb(i, nlev) = ilev - 1 ! 1
87        ELSE
88          lt(i, nlev) = 2
89          lb(i, nlev) = 1
90        END IF
91      END DO
92      DO k = 1, ilev - 1
93        DO i = 1, klon
94          pbot = pgcm(i, k)
95          ptop = pgcm(i, k+1)
96          IF (ptop<=pres(nlev) .AND. pbot>=pres(nlev)) THEN
97            lt(i, nlev) = k + 1
98            lb(i, nlev) = k
99          END IF
100        END DO
101      END DO
102
103      ! Interpolation lineaire:
104      DO i = 1, klon
105        ! interpolation en logarithme de pression:
106
107        ! ...   Modif . P. Le Van    ( 20/01/98) ....
108        ! Modif Frederic Hourdin (3/01/02)
109
110        aist(i, nlev) = log(pgcm(i,lb(i,nlev))/pres(nlev))/log(pgcm(i,lb(i, &
111          nlev))/pgcm(i,lt(i,nlev)))
112        aisb(i, nlev) = log(pres(nlev)/pgcm(i,lt(i,nlev)))/log(pgcm(i,lb(i, &
113          nlev))/pgcm(i,lt(i,nlev)))
114      END DO
115    END DO
116
117  END IF ! lnew
118
119  ! ======================================================================
120  ! inteprollation
121  ! ET je mets les vents a zero quand je rencontre une montagne
122  ! ======================================================================
123
124  DO nlev = 1, klevstd
125    DO i = 1, klon
126      IF (pgcm(i,1)<pres(nlev)) THEN
127        qpres(i, nlev) = missing_val
128      ELSE
129        qpres(i, nlev) = qgcm(i, lb(i,nlev))*aisb(i, nlev) + &
130          qgcm(i, lt(i,nlev))*aist(i, nlev)
131      END IF
132    END DO
133  END DO
134
135
136  RETURN
137END SUBROUTINE plevel_new
Note: See TracBrowser for help on using the repository browser.