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