! ! $Header$ ! c c subroutine read_pstoke(irec, . zrec,zklono,zklevo,airefi,phisfi, . t,mfu,mfd,en_u,de_u,en_d,de_d,coefh, . fm_therm,en_therm, . frac_impa,frac_nucl,pyu1,pyv1,ftsol,psrf) C****************************************************************************** C Frederic HOURDIN, Abderrahmane IDELKADI C Lecture des parametres physique stockes online necessaires pour C recalculer offline le transport de traceurs sur une grille 2x plus fine que C celle online C A FAIRE : une seule routine au lieu de 2 (lectflux, redecoupe)! C****************************************************************************** USE dimphy IMPLICIT NONE #include "netcdf.inc" #include "dimensions.h" #include "paramet.h" #include "comconst.h" #include "comgeom.h" #include "temps.h" #include "ener.h" #include "logic.h" #include "description.h" #include "serre.h" #include "indicesol.h" #include "control.h" cccc#include "dimphy.h" integer klono,klevo,imo,jmo parameter (imo=iim/2,jmo=(jjm+1)/2) parameter(klono=(jmo-1)*imo+2,klevo=llm) REAL phisfi(klono) REAL phisfi2(imo,jmo+1),airefi2(imo,jmo+1) REAL mfu(klono,klevo), mfd(klono,klevo) REAL en_u(klono,klevo), de_u(klono,klevo) REAL en_d(klono,klevo), de_d(klono,klevo) REAL coefh(klono,klevo) REAL fm_therm(klono,klevo),en_therm(klono,klevo) REAL mfu2(imo,jmo+1,klevo), mfd2(imo,jmo+1,klevo) REAL en_u2(imo,jmo+1,klevo), de_u2(imo,jmo+1,klevo) REAL en_d2(imo,jmo+1,klevo), de_d2(imo,jmo+1,klevo) REAL coefh2(imo,jmo+1,klevo) REAL fm_therm2(imo,jmo+1,klevo) REAL en_therm2(imo,jmo+1,klevo) REAL pl(klevo) integer irec integer xid,yid,zid,tid real zrec,zklono,zklevo,zim,zjm integer ncrec,ncklono,ncklevo,ncim,ncjm real airefi(klono) character*20 namedim c !! attention !! c attention il y a aussi le pb de def klono c dim de phis?? REAL frac_impa(klono,klevo), frac_nucl(klono,klevo) REAL frac_impa2(imo,jmo+1,klevo), . frac_nucl2(imo,jmo+1,klevo) REAL pyu1(klono), pyv1(klono) REAL pyu12(imo,jmo+1), pyv12(imo,jmo+1) REAL ftsol(klono,nbsrf) REAL psrf(klono,nbsrf) REAL ftsol1(klono),ftsol2(klono),ftsol3(klono),ftsol4(klono) REAL psrf1(klono),psrf2(klono),psrf3(klono),psrf4(klono) REAL ftsol12(imo,jmo+1),ftsol22(imo,jmo+1), . ftsol32(imo,jmo+1), . ftsol42(imo,jmo+1) REAL psrf12(imo,jmo+1),psrf22(imo,jmo+1),psrf32(imo,jmo+1), . psrf42(imo,jmo+1) REAL t(klono,klevo) REAL t2(imo,jmo+1,klevo) integer ncidp save ncidp integer varidt integer varidmfu, varidmfd, varidps, varidenu, variddeu integer varidend,varidded,varidch,varidfi,varidfn integer varidfmth,varidenth integer varidyu1,varidyv1,varidpl,varidai,varididvt integer varidfts1,varidfts2,varidfts3,varidfts4 integer varidpsr1,varidpsr2,varidpsr3,varidpsr4 save varidmfu, varidmfd, varidps, varidenu, variddeu save varidend,varidded,varidch,varidfi,varidfn save varidfmth,varidenth save varidyu1,varidyv1,varidpl,varidai,varididvt save varidfts1,varidfts2,varidfts3,varidfts4 save varidpsr1,varidpsr2,varidpsr3,varidpsr4 save varidt integer l, i integer start(4),count(4),status real rcode logical first save first data first/.true./ c --------------------------------------------- c Initialisation de la lecture des fichiers c --------------------------------------------- if (irec .eq. 0) then ncidp=NCOPN('phystoke.nc',NCNOWRIT,rcode) varidps=NCVID(ncidp,'phis',rcode) print*,'ncidp,varidps',ncidp,varidps varidpl=NCVID(ncidp,'sig_s',rcode) print*,'ncidp,varidpl',ncidp,varidpl varidai=NCVID(ncidp,'aire',rcode) print*,'ncidp,varidai',ncidp,varidai c A FAIRE: Es-il necessaire de stocke t? varidt=NCVID(ncidp,'t',rcode) print*,'ncidp,varidt',ncidp,varidt varidmfu=NCVID(ncidp,'mfu',rcode) print*,'ncidp,varidmfu',ncidp,varidmfu varidmfd=NCVID(ncidp,'mfd',rcode) print*,'ncidp,varidmfd',ncidp,varidmfd varidenu=NCVID(ncidp,'en_u',rcode) print*,'ncidp,varidenu',ncidp,varidenu variddeu=NCVID(ncidp,'de_u',rcode) print*,'ncidp,variddeu',ncidp,variddeu varidend=NCVID(ncidp,'en_d',rcode) print*,'ncidp,varidend',ncidp,varidend varidded=NCVID(ncidp,'de_d',rcode) print*,'ncidp,varidded',ncidp,varidded varidch=NCVID(ncidp,'coefh',rcode) print*,'ncidp,varidch',ncidp,varidch c abder (pour thermiques) varidfmth=NCVID(ncidp,'fm_th',rcode) print*,'ncidp,varidfmth',ncidp,varidfmth varidenth=NCVID(ncidp,'en_th',rcode) print*,'ncidp,varidenth',ncidp,varidenth varidfi=NCVID(ncidp,'frac_impa',rcode) print*,'ncidp,varidfi',ncidp,varidfi varidfn=NCVID(ncidp,'frac_nucl',rcode) print*,'ncidp,varidfn',ncidp,varidfn varidyu1=NCVID(ncidp,'pyu1',rcode) print*,'ncidp,varidyu1',ncidp,varidyu1 varidyv1=NCVID(ncidp,'pyv1',rcode) print*,'ncidp,varidyv1',ncidp,varidyv1 varidfts1=NCVID(ncidp,'ftsol1',rcode) print*,'ncidp,varidfts1',ncidp,varidfts1 varidfts2=NCVID(ncidp,'ftsol2',rcode) print*,'ncidp,varidfts2',ncidp,varidfts2 varidfts3=NCVID(ncidp,'ftsol3',rcode) print*,'ncidp,varidfts3',ncidp,varidfts3 varidfts4=NCVID(ncidp,'ftsol4',rcode) print*,'ncidp,varidfts4',ncidp,varidfts4 varidpsr1=NCVID(ncidp,'psrf1',rcode) print*,'ncidp,varidpsr1',ncidp,varidpsr1 varidpsr2=NCVID(ncidp,'psrf2',rcode) print*,'ncidp,varidpsr2',ncidp,varidpsr2 varidpsr3=NCVID(ncidp,'psrf3',rcode) print*,'ncidp,varidpsr3',ncidp,varidpsr3 varidpsr4=NCVID(ncidp,'psrf4',rcode) print*,'ncidp,varidpsr4',ncidp,varidpsr4 c ID pour les dimensions status = nf_inq_dimid(ncidp,'y',yid) status = nf_inq_dimid(ncidp,'x',xid) status = nf_inq_dimid(ncidp,'sig_s',zid) status = nf_inq_dimid(ncidp,'time_counter',tid) c lecture des dimensions status = nf_inq_dim(ncidp,yid,namedim,ncjm) status = nf_inq_dim(ncidp,xid,namedim,ncim) status = nf_inq_dim(ncidp,zid,namedim,ncklevo) status = nf_inq_dim(ncidp,tid,namedim,ncrec) zrec=ncrec zklevo=ncklevo zim=ncim zjm=ncjm zklono=zim*(zjm-2)+2 write(*,*) 'read_pstoke : zrec = ', zrec write(*,*) 'read_pstoke : zklevo = ', zklevo write(*,*) 'read_pstoke : zim = ', zim write(*,*) 'read_pstoke : zjm = ', zjm write(*,*) 'read_pstoke : zklono = ', zklono c niveaux de pression #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidpl,1,zklevo,pl) #else status=NF_GET_VARA_REAL(ncidp,varidpl,1,zklevo,pl) #endif c lecture de aire et phis start(1)=1 start(2)=1 start(3)=1 start(4)=0 count(1)=zim count(2)=zjm count(3)=1 count(4)=0 c phis #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidps,start,count,phisfi2) #else status=NF_GET_VARA_REAL(ncidp,varidps,start,count,phisfi2) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,phisfi2,phisfi) c aire #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidai,start,count,airefi2) #else status=NF_GET_VARA_REAL(ncidp,varidai,start,count,airefi2) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,airefi2,airefi) else print*,'ok1' c --------------------- c lecture des champs c --------------------- print*,'WARNING!!! Il n y a pas de test de coherence' print*,'sur le nombre de niveaux verticaux dans le fichier nc' start(1)=1 start(2)=1 start(3)=1 start(4)=irec count(1)=zim count(2)=zjm count(3)=zklevo count(4)=1 C *** Lessivage****************************************************** c frac_impa #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidfi,start,count,frac_impa2) #else status=NF_GET_VARA_REAL(ncidp,varidfi,start,count,frac_impa2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,frac_impa2,frac_impa) c frac_nucl #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidfn,start,count,frac_nucl2) #else status=NF_GET_VARA_REAL(ncidp,varidfn,start,count,frac_nucl2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,frac_nucl2,frac_nucl) C*** Temperature ****************************************************** c abder t #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidt,start,count,t2) #else status=NF_GET_VARA_REAL(ncidp,varidt,start,count,t2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,t2,t) C*** Flux pour le calcul de la convection TIEDTK *********************** c mfu #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidmfu,start,count,mfu2) #else status=NF_GET_VARA_REAL(ncidp,varidmfu,start,count,mfu2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,mfu2,mfu) c mfd #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidmfd,start,count,mfd2) #else status=NF_GET_VARA_REAL(ncidp,varidmfd,start,count,mfd2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,mfd2,mfd) c en_u #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidenu,start,count,en_u2) #else status=NF_GET_VARA_REAL(ncidp,varidenu,start,count,en_u2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,en_u2,en_u) c de_u #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,variddeu,start,count,de_u2) #else status=NF_GET_VARA_REAL(ncidp,variddeu,start,count,de_u2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,de_u2,de_u) c en_d #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidend,start,count,en_d2) #else status=NF_GET_VARA_REAL(ncidp,varidend,start,count,en_d2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,en_d2,en_d) c de_d #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidded,start,count,de_d2) #else status=NF_GET_VARA_REAL(ncidp,varidded,start,count,de_d2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,de_d2,de_d) C **** Coeffecient du mellange turbulent********************************** c coefh #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidch,start,count,coefh2) #else status=NF_GET_VARA_REAL(ncidp,varidch,start,count,coefh2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,coefh2,coefh) C*** Flux ascendant et entrant pour les Thermiques************************ cabder thermiques #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidfmth,start,count,fm_therm2) #else status=NF_GET_VARA_REAL(ncidp,varidfmth,start,count,fm_therm2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,fm_therm2,fm_therm) #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidenth,start,count,en_therm2) #else status=NF_GET_VARA_REAL(ncidp,varidenth,start,count,en_therm2) #endif call gr_ecrit_fi(klevo,klono,imo,jmo+1,en_therm2,en_therm) C*** Vitesses aux sol ****************************************************** start(3)=irec start(4)=0 count(3)=1 count(4)=0 c pyu1 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidyu1,start,count,pyu12) #else status=NF_GET_VARA_REAL(ncidp,varidyu1,start,count,pyu12) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,pyu12,pyu1) c pyv1 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidyv1,start,count,pyv12) #else status=NF_GET_VARA_REAL(ncidp,varidyv1,start,count,pyv12) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,pyv12,pyv1) C*** Temperature au sol ******************************************** c ftsol1 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidfts1,start,count,ftsol12) #else status=NF_GET_VARA_REAL(ncidp,varidfts1,start,count,ftsol12) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,ftsol12,ftsol1) c ftsol2 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidfts2,start,count,ftsol22) #else status=NF_GET_VARA_REAL(ncidp,varidfts2,start,count,ftsol22) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,ftsol22,ftsol2) c ftsol3 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidfts3,start,count,ftsol32) #else status=NF_GET_VARA_REAL(ncidp,varidfts3,start,count,ftsol32) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,ftsol32,ftsol3) c ftsol4 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidfts4,start,count,ftsol42) #else status=NF_GET_VARA_REAL(ncidp,varidfts4,start,count,ftsol42) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,ftsol42,ftsol4) C*** Nature du sol ************************************************** c psrf1 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidpsr1,start,count,psrf12) #else status=NF_GET_VARA_REAL(ncidp,varidpsr1,start,count,psrf12) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,psrf12,psrf1) c psrf2 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidpsr2,start,count,psrf22) #else status=NF_GET_VARA_REAL(ncidp,varidpsr2,start,count,psrf22) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,psrf22,psrf2) c psrf3 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidpsr3,start,count,psrf32) #else status=NF_GET_VARA_REAL(ncidp,varidpsr3,start,count,psrf32) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,psrf32,psrf3) c psrf4 #ifdef NC_DOUBLE status=NF_GET_VARA_DOUBLE(ncidp,varidpsr4,start,count,psrf42) #else status=NF_GET_VARA_REAL(ncidp,varidpsr4,start,count,psrf42) #endif call gr_ecrit_fi(1,klono,imo,jmo+1,psrf42,psrf4) do i = 1,klono psrf(i,1) = psrf1(i) psrf(i,2) = psrf2(i) psrf(i,3) = psrf3(i) psrf(i,4) = psrf4(i) ftsol(i,1) = ftsol1(i) ftsol(i,2) = ftsol2(i) ftsol(i,3) = ftsol3(i) ftsol(i,4) = ftsol4(i) enddo endif return end