Index: trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/readhead_NC.F =================================================================== --- trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/readhead_NC.F (revision 1974) +++ trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/readhead_NC.F (revision 1977) @@ -5,5 +5,5 @@ USE comvert_mod, ONLY: aps,bps,preff USE comconst_mod, ONLY: im,jm,lllm,daysec,dtvr, - . rad,omeg,g,cpp,kappa,r + . rad,omeg,g,cpp,kappa,r,pi USE temps_mod, ONLY: day_ini USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0 @@ -110,7 +110,7 @@ ENDIF - ierr = NF_INQ_VARID (nid, "rlonu", nvarid) - IF (ierr .NE. NF_NOERR) THEN - PRINT*, "readhead_NC: Le champ est absent" + ierr = NF_INQ_VARID (nid, "longitude", nvarid) + IF (ierr .NE. NF_NOERR) THEN + PRINT*, "readhead_NC: Le champ est absent" CALL abort ENDIF @@ -121,11 +121,13 @@ #endif IF (ierr .NE. NF_NOERR) THEN - PRINT*, "readhead_NC: Lecture echouee pour " - CALL abort - ENDIF - - ierr = NF_INQ_VARID (nid, "rlatv", nvarid) - IF (ierr .NE. NF_NOERR) THEN - PRINT*, "readhead_NC: Le champ est absent" + PRINT*, "readhead_NC: Lecture echouee pour " + CALL abort + ENDIF +! convert it to radians + rlonu(:)=rlonu(:)*pi/180. + + ierr = NF_INQ_VARID (nid, "latitude", nvarid) + IF (ierr .NE. NF_NOERR) THEN + PRINT*, "readhead_NC: Le champ est absent" CALL abort ENDIF @@ -136,13 +138,15 @@ #endif IF (ierr .NE. NF_NOERR) THEN - PRINT*, "readhead_NC: Lecture echouee pour rlatv" - CALL abort - ENDIF - - ierr = NF_GET_VAR_REAL(nid, nvarid, cv) - IF (ierr .NE. NF_NOERR) THEN - PRINT*, "readhead_NC: Lecture echouee pour " - CALL abort - ENDIF + PRINT*, "readhead_NC: Lecture echouee pour latitude" + CALL abort + ENDIF +! convert it to radians + rlatv(:)=rlatv(:)*pi/180. + +! ierr = NF_GET_VAR_REAL(nid, nvarid, cv) +! IF (ierr .NE. NF_NOERR) THEN +! PRINT*, "readhead_NC: Lecture echouee pour " +! CALL abort +! ENDIF c c Lecture des aires des mailles: Index: trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/xvik.F =================================================================== --- trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/xvik.F (revision 1974) +++ trunk/LMDZ.MARS/libf/dynphy_lonlat/phymars/xvik.F (revision 1977) @@ -15,10 +15,10 @@ -#include "dimensions.h" -#include "paramet.h" -#include "comdissip.h" -#include "comgeom2.h" + include "dimensions.h" + include "paramet.h" + include "comdissip.h" + include "comgeom2.h" !#include "control.h" -#include "netcdf.inc" + include "netcdf.inc" @@ -86,5 +86,6 @@ c ---------------- - unanj=667.9 + chr2="0" + unanj=669. print*,'WARNING!!!',unanj,'Jours/an' nc=.true. @@ -96,5 +97,5 @@ c lecture du fichier xvik.def - phivik(1)=-3627 + phivik(1)=-3637 phivik(2)=-4505 @@ -179,7 +180,8 @@ WRITE (*,*) 'day0 = ' , day0 + CALL conf_gcm( 99, .TRUE. ) CALL iniconst CALL inigeom - CALL inifilr +! CALL inifilr @@ -205,5 +207,5 @@ ENDIF - PRINT*,'temps',(temps(itau),itau=1,10) + PRINT*,'temps(1:10)',(temps(itau),itau=1,10) c----------------------------------------------------------------------- @@ -212,5 +214,5 @@ latvik(1)=22.27*pi/180. - lonvik(1)=-47.9*pi/180. + lonvik(1)=-47.95*pi/180. latvik(2)=47.67*pi/180. lonvik(2)=(360.-225.71)*pi/180. @@ -218,15 +220,25 @@ c ponderations pour les 4 points autour de Viking DO iv=1,2 - CALL coordij(lonvik(iv),latvik(iv),ivik(iv),jvik(iv)) - IF(lonvik(iv).lt.rlonv(ivik(iv))) THEN - ivik(iv)=ivik(iv)-1 - ENDIF - IF(latvik(iv).gt.rlatu(jvik(iv))) THEN - jvik(iv)=jvik(iv)-1 - ENDIF - zalpha=(lonvik(iv)-rlonv(ivik(iv)))/ - s (rlonv(ivik(iv)+1)-rlonv(ivik(iv))) - zbeta=(latvik(iv)-rlatu(jvik(iv)))/ - s (rlatu(jvik(iv)+1)-rlatu(jvik(iv))) + ! locate index of GCM grid points near VL + do i=1,iim + ! we know longitudes are ordered -180...180 + if ((lonvik(iv).ge.rlonu(i)).and. + & (lonvik(iv).le.rlonu(i+1))) then + ivik(iv)=i + exit + endif + enddo + do j=1,jjm-1 + !we know tha latitudes are ordered 90...-90 + if ((latvik(iv).le.rlatv(j)).and. + & (latvik(iv).ge.rlatv(j+1))) then + jvik(iv)=j + exit + endif + enddo + zalpha=(lonvik(iv)-rlonu(ivik(iv)))/ + s (rlonu(ivik(iv)+1)-rlonu(ivik(iv))) + zbeta=(latvik(iv)-rlatv(jvik(iv)))/ + s (rlatv(jvik(iv)+1)-rlatv(jvik(iv))) zw(0,0,iv)=(1.-zalpha)*(1.-zbeta) zw(1,0,iv)=zalpha*(1.-zbeta) @@ -346,5 +358,5 @@ PRINT*,'day ',day CALL solarlong(day+day0,sollong) - sol=day+day0+461. + sol=day+day0+461. ! aded offset to sync with VL mission "sol 1" iyear=sol/unanj WRITE (*,*) 'iyear',iyear @@ -396,5 +408,5 @@ OPEN(97,file='xprestot'//chr2,form='formatted') c OPEN(98,file='xlat37_'//chr2,form='formatted') - WRITE(98,'(f5.1,16f7.1)') 0.,(rlonv(i)*180./pi,i=1,iim,4) +c WRITE(98,'(f5.1,16f7.1)') 0.,(rlonv(i)*180./pi,i=1,iim,4) ENDIF @@ -433,6 +445,6 @@ IF(.NOT.firstcal) THEN - WRITE(98,'(f5.1,16f7.3)') - s dayw,(ps(i,37),i=1,iim,4) +c WRITE(98,'(f5.1,16f7.3)') +c s dayw,(ps(i,37),i=1,iim,4) c boucle sur les sites vikings: @@ -448,6 +460,6 @@ zp2_sm=0. zt=0. - zu=0. - zv=0. +! zu=0. +! zv=0. DO jj=0,1 j=jvik(iv)+jj @@ -459,6 +471,6 @@ zp1=zp1+zw(ii,jj,iv)*log(ps(i,j)) ! interpolate in log(P) WRITE (*,*) 'ps autour iv',ps(i,j),iv - zu=zu+zw(ii,jj,iv)*unat(i,j,1)/cu(i,j) - zv=zv+zw(ii,jj,iv)*vnat(i,j,1)/cv(i,j) +! zu=zu+zw(ii,jj,iv)*unat(i,j,1)/cu(i,j) +! zv=zv+zw(ii,jj,iv)*vnat(i,j,1)/cv(i,j) ENDDO ENDDO @@ -481,10 +493,10 @@ c sorties eventuelles de vent - IF(lvent) THEN - WRITE(ifile(iv)+16,'(2e15.5)') - s dayw,zu - WRITE(ifile(iv)+12,'(2e15.5)') - s dayw,zv - ENDIF +! IF(lvent) THEN +! WRITE(ifile(iv)+16,'(2e15.5)') +! s dayw,zu +! WRITE(ifile(iv)+12,'(2e15.5)') +! s dayw,zv +! ENDIF ENDDO c IF (lcal) THEN