| 1 | SUBROUTINE load_psi( |
|---|
| 2 | S psurf, ztop, ksive, |
|---|
| 3 | S temp, psimap, deltapsimap) |
|---|
| 4 | |
|---|
| 5 | use dimphy |
|---|
| 6 | IMPLICIT none |
|---|
| 7 | |
|---|
| 8 | #include "dimensions.h" |
|---|
| 9 | #include "YOMCST.h" |
|---|
| 10 | #include "comcstVE.h" |
|---|
| 11 | C |
|---|
| 12 | C ------------------------------------------------------------------ |
|---|
| 13 | C |
|---|
| 14 | C PURPOSE. |
|---|
| 15 | C -------- |
|---|
| 16 | C |
|---|
| 17 | c This routine loads the longwave matrix of factors Ksi |
|---|
| 18 | c (interpolated for a given column) |
|---|
| 19 | c to build the Net Exchange Rates matrix Psi. |
|---|
| 20 | c Psi(i,j,nu) = Ksi(i,j,nu) * ( B(i,nu)-B(j,nu) ) |
|---|
| 21 | c |
|---|
| 22 | c The Ksi matrixes have been computed by Vincent Eymet |
|---|
| 23 | C |
|---|
| 24 | c The NER matrix is then integrated in frequency. |
|---|
| 25 | c |
|---|
| 26 | C AUTHOR. |
|---|
| 27 | C ------- |
|---|
| 28 | C Sebastien Lebonnois |
|---|
| 29 | C |
|---|
| 30 | C MODIFICATIONS. |
|---|
| 31 | C -------------- |
|---|
| 32 | C version multimatrice (topographie, sommet nuages): 20/12/2006 |
|---|
| 33 | C ------------------------------------------------------------------ |
|---|
| 34 | C |
|---|
| 35 | C* ARGUMENTS: |
|---|
| 36 | C |
|---|
| 37 | c inputs |
|---|
| 38 | real psurf(klon) ! Surface pressure |
|---|
| 39 | real ztop(klon) ! Altitude of the top of cloud deck (km) |
|---|
| 40 | real ksive(0:kflev+1,0:kflev+1,nnuve,nbmat) ! ksi matrixes in Vincent's file |
|---|
| 41 | real temp(klon,0:kflev+1) ! Temperature in layer (K) |
|---|
| 42 | c outputs |
|---|
| 43 | real psimap(0:kflev+1,0:kflev+1,klon) |
|---|
| 44 | real deltapsimap(0:kflev+1,0:kflev+1,klon) |
|---|
| 45 | |
|---|
| 46 | c local variables |
|---|
| 47 | integer i,j,ig,band,nlve |
|---|
| 48 | integer mat,m,mat0 |
|---|
| 49 | character*100 file |
|---|
| 50 | real bplck(0:kflev+1,nnuve) ! Planck luminances in table layers |
|---|
| 51 | real y(0:kflev,nnuve) ! intermediaire Planck |
|---|
| 52 | real zdblay(0:kflev+1,nnuve) ! gradient en temperature de planck |
|---|
| 53 | real ksi |
|---|
| 54 | real factp,factz |
|---|
| 55 | |
|---|
| 56 | nlve = kflev ! (doit correspondre, pour bplck, y et zdblay) |
|---|
| 57 | |
|---|
| 58 | c ----------------------- |
|---|
| 59 | c Main loop on grid point |
|---|
| 60 | c ----------------------- |
|---|
| 61 | do ig=1,klon |
|---|
| 62 | |
|---|
| 63 | c Planck function |
|---|
| 64 | c --------------- |
|---|
| 65 | |
|---|
| 66 | do band=1,nnuve |
|---|
| 67 | do j=0,nlve |
|---|
| 68 | c B(T,l) = al/(exp(bl/T)-1) |
|---|
| 69 | y(j,band) = exp(bl(band)/temp(ig,j))-1. |
|---|
| 70 | bplck(j,band) = al(band)/(y(j,band)) |
|---|
| 71 | zdblay(j,band) = al(band)*bl(band)*exp(bl(band)/temp(ig,j))/ |
|---|
| 72 | . ((temp(ig,j)*temp(ig,j))*(y(j,band)*y(j,band))) |
|---|
| 73 | enddo |
|---|
| 74 | bplck(nlve+1,band) = 0.0 |
|---|
| 75 | zdblay(nlve+1,band)= 0.0 |
|---|
| 76 | enddo |
|---|
| 77 | |
|---|
| 78 | c interpolating ksi |
|---|
| 79 | c and |
|---|
| 80 | c computing psi and deltapsi |
|---|
| 81 | c --------------------------- |
|---|
| 82 | |
|---|
| 83 | c init |
|---|
| 84 | do j=0,nlve+1 |
|---|
| 85 | do i=0,nlve+1 |
|---|
| 86 | psimap(i,j,ig) = 0.0 ! positif quand nrj de i->j |
|---|
| 87 | deltapsimap(i,j,ig) = 0.0 |
|---|
| 88 | enddo |
|---|
| 89 | enddo |
|---|
| 90 | |
|---|
| 91 | c finding the right matrixes |
|---|
| 92 | mat0 = 0 |
|---|
| 93 | do mat=1,nbmat-nbztopve |
|---|
| 94 | if ( (psurfve(mat).ge.psurf(ig)) |
|---|
| 95 | . .and.(psurfve(mat+nbztopve).lt.psurf(ig)) |
|---|
| 96 | . .and.(ztopve(mat).lt.ztop(ig)) |
|---|
| 97 | . .and.(ztopve(mat+1).ge.ztop(ig)) ) then |
|---|
| 98 | mat0 = mat |
|---|
| 99 | c print*,'ig=',ig,' mat0=',mat |
|---|
| 100 | factp = (psurf(ig) -psurfve(mat)) |
|---|
| 101 | . /(psurfve(mat+nbztopve)-psurfve(mat)) |
|---|
| 102 | factz = (ztop(ig) -ztopve(mat)) |
|---|
| 103 | . /(ztopve(mat+1)-ztopve(mat)) |
|---|
| 104 | endif |
|---|
| 105 | enddo |
|---|
| 106 | if (mat0.eq.0) then |
|---|
| 107 | write(*,*) 'This is subroutine load_psi' |
|---|
| 108 | print*,'Probleme pour interpolation au point ig=',ig |
|---|
| 109 | print*,'psurf = ',psurf(ig),' ztop = ',ztop(ig) |
|---|
| 110 | stop |
|---|
| 111 | endif |
|---|
| 112 | |
|---|
| 113 | c------------------ |
|---|
| 114 | c !!TEST!! Matrice unique fixee: psurf = 90 bar, ztop = 70 |
|---|
| 115 | c mat0 = 24 |
|---|
| 116 | c print*,'MATRICE UNIQUE: ',mat0,' / ps=',psurfve(mat0), |
|---|
| 117 | c . ' / ztop=',ztopve(mat0) |
|---|
| 118 | c------------------ |
|---|
| 119 | |
|---|
| 120 | c interpolation of ksi and computation of psi,deltapsi |
|---|
| 121 | do band=1,nnuve |
|---|
| 122 | do j=0,nlve+1 |
|---|
| 123 | do i=0,nlve+1 |
|---|
| 124 | ksi = ksive(i,j,band,mat0)*(1-factz)*(1-factp) |
|---|
| 125 | . +ksive(i,j,band,mat0+1)*factz *(1-factp) |
|---|
| 126 | . +ksive(i,j,band,mat0+nbztopve)*(1-factz)*factp |
|---|
| 127 | . +ksive(i,j,band,mat0+nbztopve+1)*factz *factp |
|---|
| 128 | psimap(i,j,ig) = psimap(i,j,ig) + |
|---|
| 129 | . ksi*(bplck(i,band)-bplck(j,band)) |
|---|
| 130 | deltapsimap(i,j,ig) = deltapsimap(i,j,ig) + |
|---|
| 131 | . ksi*zdblay(i,band) |
|---|
| 132 | enddo |
|---|
| 133 | enddo |
|---|
| 134 | enddo |
|---|
| 135 | |
|---|
| 136 | enddo !ig |
|---|
| 137 | c ----------------------- |
|---|
| 138 | c End loop on grid point |
|---|
| 139 | c ----------------------- |
|---|
| 140 | |
|---|
| 141 | c print*,"LOAD_PSI OK" |
|---|
| 142 | |
|---|
| 143 | return |
|---|
| 144 | end |
|---|
| 145 | |
|---|
