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