Changeset 2908
- Timestamp:
- 2011-10-13T07:52:23+02:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_r2802_MERCATOR10_diadct/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90
r2864 r2908 21 21 !!---------------------------------------------------------------------- 22 22 !!---------------------------------------------------------------------- 23 !! dia_dct : compute the transport t rough a sec.23 !! dia_dct : compute the transport through a sec. 24 24 !! dia_dct_init : read namelist. 25 25 !! readsec : read sections description and pathway … … 78 78 79 79 TYPE SECTION 80 CHARACTER(len=60) :: name ! name of the sec81 LOGICAL :: llstrpond ! true if you want the computation of salinityand82 ! temperature balanced by the transport83 LOGICAL :: ll_ice_section ! icesurf and icevol computation84 LOGICAL :: ll_date_line ! = T if the section crosses the date-line85 TYPE(COORD_SECTION), DIMENSION(2) :: coordSec ! longitude and latitude of the extremities of the sec86 INTEGER :: nb_class ! number of boundaries for density classes87 INTEGER, DIMENSION(nb_point_max) :: direction ! vector direction of the point in the section88 CHARACTER(len=40),DIMENSION(nb_class_max) :: classname ! caracteristics of the class89 REAL(wp), DIMENSION(nb_class_max) :: zsigi ,&! in-situ density classes (99 if you don't want)90 zsigp ,&! potential density classes (99 if you don't want)91 zsal ,&! salinity classes (99 if you don't want)92 ztem ,&! temperature classes(99 if you don't want)93 zlay ! level classes (99 if you don't want)80 CHARACTER(len=60) :: name ! name of the sec 81 LOGICAL :: llstrpond ! true if you want the computation of salt and 82 ! heat transports 83 LOGICAL :: ll_ice_section ! icesurf and icevol computation 84 LOGICAL :: ll_date_line ! = T if the section crosses the date-line 85 TYPE(COORD_SECTION), DIMENSION(2) :: coordSec ! longitude and latitude of the extremities of the sec 86 INTEGER :: nb_class ! number of boundaries for density classes 87 INTEGER, DIMENSION(nb_point_max) :: direction ! vector direction of the point in the section 88 CHARACTER(len=40),DIMENSION(nb_class_max) :: classname ! caracteristics of the class 89 REAL(wp), DIMENSION(nb_class_max) :: zsigi ,&! in-situ density classes (99 if you don't want) 90 zsigp ,&! potential density classes (99 if you don't want) 91 zsal ,&! salinity classes (99 if you don't want) 92 ztem ,&! temperature classes(99 if you don't want) 93 zlay ! level classes (99 if you don't want) 94 94 REAL(wp), DIMENSION(nb_type_class,nb_class_max) :: transport ! transport output 95 REAL(wp) :: slopeSection ! coeff directeur de la section96 INTEGER :: nb_point ! nb de points de la section97 TYPE(POINT_SECTION),DIMENSION(nb_point_max) :: listPoint ! list of point in the section95 REAL(wp) :: slopeSection ! section's slopesection 96 INTEGER :: nb_point ! section's number of points 97 TYPE(POINT_SECTION),DIMENSION(nb_point_max) :: listPoint ! section's list point 98 98 END TYPE SECTION 99 99 … … 161 161 LOGICAL :: lldebug =.FALSE. !debug a section 162 162 CHARACTER(len=160) :: clfileout !fileout name 163 164 165 INTEGER , DIMENSION(1):: ish ! tmp array for mpp_sum 166 INTEGER , DIMENSION(3):: ish2 ! " 167 REAL(wp), DIMENSION(nb_sec_max*nb_type_class*nb_class_max):: zwork ! " 168 REAL(wp), DIMENSION(nb_sec_max,nb_type_class,nb_class_max):: zsum ! " 169 163 170 !!--------------------------------------------------------------------- 164 171 … … 183 190 !Compute transport through section 184 191 CALL transport(secs(jsec),lldebug) 192 193 ENDDO 185 194 186 187 188 195 IF( MOD(kt,nn_dctwri)==0 )THEN 196 197 IF( lwp .AND. kt==nit000+nn_dctwri-1 )WRITE(numout,*)" diadct: write at kt = ",kt 189 198 190 !Write the transport 199 !Sum on all procs 200 IF( lk_mpp )THEN 201 ish(1) = nb_sec_max*nb_type_class*nb_class_max 202 ish2 = (/nb_sec_max,nb_type_class,nb_class_max/) 203 DO jsec=1,nb_sec ; zsum(jsec,:,:) = secs(jsec)%transport(:,:) ; ENDDO 204 zwork(:)= RESHAPE(zsum(:,:,:), ish ) 205 CALL mpp_sum(zwork, ish(1)) 206 zsum(:,:,:)= RESHAPE(zwork,ish2) 207 DO jsec=1,nb_sec ; secs(jsec)%transport(:,:) = zsum(jsec,:,:) ; ENDDO 208 ENDIF 209 210 !Write the transport 211 DO jsec=1,nb_sec 212 191 213 IF( lwp )CALL dia_dct_wri(kt,jsec,secs(jsec)) 192 214 … … 194 216 secs(jsec)%transport(:,:)=0. 195 217 196 ENDIF 197 ENDDO 218 ENDDO 219 220 ENDIF 221 198 222 ENDIF 199 223 … … 310 334 ENDDO 311 335 ENDIF 312 !?IF( jsec==nn_secdebug .OR. nn_secdebug==-1 )THEN313 ! DO jpt=1,iptglo314 ! WRITE(numout,*)' # I J ',jpt,coordtemp(jpt)315 ! ENDDO316 !ENDIF317 336 318 337 !Now each proc selects only points that are in its domain: … … 335 354 secs(jsec)%listPoint(iptloc) = POINT_SECTION(mi0(iiglo),mj0(ijglo)) ! store local coordinates 336 355 secs(jsec)%direction(iptloc) = directemp(jpt) ! store local direction 337 !WRITE(narea+200,*)' # I J : ',iiloc,ijloc338 356 ENDIF 339 357 … … 356 374 357 375 !remove redundant points between processors 358 !look NEMO documentation for more explanations 359 ! lldebug = .FALSE. ; IF ( (jsec==nn_secdebug .OR. nn_secdebug==-1) .AND. lwp ) lldebug = .TRUE. 360 lldebug = .FALSE. ; IF ( (jsec==nn_secdebug .OR. nn_secdebug==-1) ) lldebug = .TRUE. 376 !------------------------------------------ 377 lldebug = .FALSE. ; IF ( (jsec==nn_secdebug .OR. nn_secdebug==-1) .AND. lwp ) lldebug = .TRUE. 361 378 IF( iptloc .NE. 0 )THEN 362 379 CALL removepoints(secs(jsec),'I','top_list',lldebug) … … 409 426 !where points could be redundant 410 427 isgn ,& !way of course in listpoint 428 ipoint ,& !way of course in listpoint 411 429 istart,iend !first and last points selected in listpoint 412 430 INTEGER :: jpoint =0 !loop on list points … … 443 461 444 462 jpoint=iextr+isgn 445 DO WHILE ( jpoint .GE. 1 .AND. jpoint .LE. sec%nb_point.AND. &446 icoord( iind,jpoint-isgn ) == itest .AND. icoord( iind,jpoint ) == itest)447 463 DO WHILE( jpoint .GE. 1 .AND. jpoint .LE. sec%nb_point .AND. & 464 icoord( iind,jpoint-isgn ) == itest .AND. icoord( iind,jpoint ) == itest ) 465 jpoint=jpoint+isgn 448 466 ENDDO 467 449 468 IF( cdextr=='bot_list')THEN ; istart=jpoint-1 ; iend=sec%nb_point 450 469 ELSE ; istart=1 ; iend=jpoint+1 … … 489 508 !! * Local variables 490 509 INTEGER :: jk,jseg,jclass, &!loop on level/segment/classes 491 iptegrid, &!=0 if section is constant along i/j492 510 isgnu , isgnv ! 493 511 INTEGER :: ii, ij ! local integer … … 506 524 507 525 TYPE(POINT_SECTION) :: k 508 INTEGER,DIMENSION(1):: ish509 INTEGER,DIMENSION(2):: ish2510 526 REAL(wp),DIMENSION(nb_type_class,nb_class_max)::zsum 511 REAL(wp),DIMENSION(nb_type_class*nb_class_max)::zwork ! temporary vector for mpp_sum512 527 !!-------------------------------------------------------- 513 528 … … 520 535 zice_surf_neg = 0._wp ; zice_surf_pos = 0._wp 521 536 zice_vol_pos = 0._wp ; zice_vol_neg = 0._wp 522 523 !compute iptegrid: iptegrid=0 if section is constant along i or j524 iptegrid = 0525 jseg = 1526 DO WHILE( (jseg .LT. MAX(sec%nb_point-1,0)) .AND. iptegrid .NE. 1 )527 IF(sec%direction(jseg) .NE. sec%direction(1)) iptegrid = 1528 jseg = jseg+1529 ENDDO530 IF(lk_mpp) CALL mpp_sum(iptegrid)531 537 532 538 !---------------------------! … … 623 629 !TEST ON THE DENSITY/SALINITY/TEMPERATURE/LEVEL! 624 630 !----------------------------------------------! 625 IF(ld_debug)write(narea+200,*)k ,jk626 631 627 632 IF ( ( ((( zrhop .GE. (sec%zsigp(jclass)+1000. )) .AND. & … … 663 668 ENDIF 664 669 #endif 665 !IF(ld_debug)write(narea+200,*)k ,jk,zTnorm666 670 !COMPUTE TRANSPORT 667 671 !zTnorm=transport through one cell for one class … … 771 775 ENDIF !end of nb_inmesh=0 case 772 776 773 !------------------!774 !SUM ON ALL PROCS !775 !------------------!776 IF( lk_mpp )THEN777 ish(1) = nb_type_class*nb_class_max ; ish2(1)=nb_type_class ; ish2(2)=nb_class_max778 zwork(:)= RESHAPE(zsum(:,:), ish )779 CALL mpp_sum(zwork, ish(1))780 zsum(:,:)= RESHAPE(zwork,ish2)781 ENDIF782 783 777 !-------------------------------| 784 778 !FINISH COMPUTING TRANSPORTS | … … 838 832 !! Method: 839 833 !! 1. Write volume transports in "volume_transport" 840 !! Unit: Sv : Surface* Velocity / 1.e6834 !! Unit: Sv : area * Velocity / 1.e6 841 835 !! 842 836 !! 2. Write heat transports in "heat_transport" 843 !! Unit: Peta W : Surface* Velocity * T * rhau * Cp / 1.e15837 !! Unit: Peta W : area * Velocity * T * rhau * Cp / 1.e15 844 838 !! 845 839 !! 3. Write salt transports in "salt_transport" 846 !! Unit: Mega m^3 / s : Surface* Velocity * S / 1.e6840 !! Unit: 10^9 g m^3 / s : area * Velocity * S / 1.e6 847 841 !! 848 842 !!------------------------------------------------------------- 849 843 !!arguments 850 INTEGER, INTENT(IN) :: kt 851 TYPE(SECTION), INTENT(INOUT) :: sec 852 INTEGER ,INTENT(IN) :: jsec 844 INTEGER, INTENT(IN) :: kt ! time-step 845 TYPE(SECTION), INTENT(INOUT) :: sec ! section to write 846 INTEGER ,INTENT(IN) :: jsec ! section number 853 847 854 848 !!local declarations 855 849 REAL(wp) ,DIMENSION(nb_type_class):: zsumclass 856 INTEGER :: jcl,ji !Dummy loop 857 CHARACTER(len=2) :: classe 858 REAL(wp) :: zbnd1,zbnd2 850 INTEGER :: jcl,ji ! Dummy loop 851 CHARACTER(len=2) :: classe ! Classname 852 REAL(wp) :: zbnd1,zbnd2 ! Class bounds 853 REAL(wp) :: zslope ! section's slope coeff 859 854 !!------------------------------------------------------------- 860 855 861 856 zsumclass(:)=0._wp 862 857 zslope = sec%slopeSection 858 859 863 860 DO jcl=1,MAX(1,sec%nb_class-1) 864 861 … … 900 897 ENDIF 901 898 !temperature classes transports 902 IF( ( sec%ztem(jcl+1) .NE. 99. ) .AND. &903 ( sec%ztem(jcl+1) .NE. 99. ) ) THEN899 IF( ( sec%ztem(jcl+1) .NE. 99._wp ) .AND. & 900 ( sec%ztem(jcl+1) .NE. 99._wp ) ) THEN 904 901 classe = 'T ' 905 902 zbnd1 = sec%ztem(jcl) … … 908 905 909 906 !write volume transport per class 910 WRITE(numdct_vol,118) ndastp,kt,jsec,sec%name, sec%slopeSection, &907 WRITE(numdct_vol,118) ndastp,kt,jsec,sec%name,zslope, & 911 908 jcl,classe,zbnd1,zbnd2,& 912 909 sec%transport(1,jcl),sec%transport(2,jcl), & … … 916 913 917 914 !write heat transport per class: 918 WRITE(numdct_heat,119) ndastp,kt,jsec,sec%name, sec%slopeSection, &915 WRITE(numdct_heat,119) ndastp,kt,jsec,sec%name,zslope, & 919 916 jcl,classe,zbnd1,zbnd2,& 920 917 sec%transport(7,jcl)*1000._wp*rcp/1.e15,sec%transport(8,jcl)*1000._wp*rcp/1.e15, & 921 918 ( sec%transport(7,jcl)+sec%transport(8,jcl) )*1000._wp*rcp/1.e15 922 919 !write salt transport per class 923 WRITE(numdct_salt,119) ndastp,kt,jsec,sec%name, sec%slopeSection, &920 WRITE(numdct_salt,119) ndastp,kt,jsec,sec%name,zslope, & 924 921 jcl,classe,zbnd1,zbnd2,& 925 922 sec%transport(9,jcl)*1000._wp/1.e9,sec%transport(10,jcl)*1000._wp/1.e9,& … … 934 931 935 932 !write total volume transport 936 WRITE(numdct_vol,118) ndastp,kt,jsec,sec%name, sec%slopeSection, &933 WRITE(numdct_vol,118) ndastp,kt,jsec,sec%name,zslope, & 937 934 jcl,"total",zbnd1,zbnd2,& 938 935 zsumclass(1),zsumclass(2),zsumclass(1)+zsumclass(2) … … 941 938 942 939 !write total heat transport 943 WRITE(numdct_heat,119) ndastp,kt,jsec,sec%name, sec%slopeSection, &940 WRITE(numdct_heat,119) ndastp,kt,jsec,sec%name,zslope, & 944 941 jcl,"total",zbnd1,zbnd2,& 945 zsumclass(7)* 1000._wp*rcp/1.e15,zsumclass(8)* 1000._wp*rcp/1.e15,(zsumclass(7)+zsumclass(8) )* 1000._wp*rcp/1.e15 942 zsumclass(7)* 1000._wp*rcp/1.e15,zsumclass(8)* 1000._wp*rcp/1.e15,& 943 (zsumclass(7)+zsumclass(8) )* 1000._wp*rcp/1.e15 946 944 !write total salt transport 947 WRITE(numdct_salt,119) ndastp,kt,jsec,sec%name, sec%slopeSection, &945 WRITE(numdct_salt,119) ndastp,kt,jsec,sec%name,zslope, & 948 946 jcl,"total",zbnd1,zbnd2,& 949 zsumclass(9)*1000._wp/1.e9,zsumclass(10)*1000._wp/1.e9,(zsumclass(9)+zsumclass(10))*1000._wp/1.e9 947 zsumclass(9)*1000._wp/1.e9,zsumclass(10)*1000._wp/1.e9,& 948 (zsumclass(9)+zsumclass(10))*1000._wp/1.e9 950 949 ENDIF 951 950 … … 953 952 IF ( sec%ll_ice_section) THEN 954 953 !write total ice volume transport 955 WRITE(numdct_vol,118) ndastp,kt,jsec,sec%name, sec%slopeSection,&954 WRITE(numdct_vol,118) ndastp,kt,jsec,sec%name,zslope,& 956 955 jcl,"ice_vol",zbnd1,zbnd2,& 957 956 sec%transport(9,1),sec%transport(10,1),& 958 957 sec%transport(9,1)+sec%transport(10,1) 959 958 !write total ice surface transport 960 WRITE(numdct_vol,118) ndastp,kt,jsec,sec%name, sec%slopeSection,&959 WRITE(numdct_vol,118) ndastp,kt,jsec,sec%name,zslope,& 961 960 jcl,"ice_surf",zbnd1,zbnd2,& 962 961 sec%transport(11,1),sec%transport(12,1), & … … 972 971 !!---------------------------------------------------------------------- 973 972 !! 974 !! Purpose: compute Temperature/Salinity/density on U-point or V-point 973 !! Purpose: compute Temperature/Salinity/density at U-point or V-point 974 !! -------- 975 975 !! 976 !! | I | I+1 | 977 !! | | | 978 !! -------------------------------- 1.Compute Tbis: 979 !! | | | interpolate T(I,J,K) et T(I,J,K+1) 980 !! | | | 981 !! | | | 2. 1.Compute Temp/Salt/Rho à U point: 982 !! K | T | | interpolate Tbis et T(I+U,J,K+1) 983 !! | | | 984 !! | | | 985 !! | | | 986 !! -------------------------------- 987 !! | | | 988 !! | | | 989 !! | | | 990 !!K+1 | Tbis U T | 991 !! | | | 992 !! | T | | 993 !! | |------------| 994 !! | | partials | 995 !! | | steps | 996 !! -------------------------------- 976 !! Method: 977 !! ------ 978 !! 979 !! ====> full step and partial step 980 !! 981 !! 982 !! | I | I+1 | Z=Temperature/Salinity/density at U-poinT 983 !! | | | 984 !! ---------------------------------------- 1. Veritcale interpolation: compute zbis 985 !! | | | interpolation between ptab(I,J,K) and ptab(I,J,K+1) 986 !! | | | zbis = 987 !! | | | [ e3w(I+1,J,K)*ptab(I,J,K) + ( e3w(I,J,K) - e3w(I+1,J,K) ) * ptab(I,J,K-1) ] 988 !! | | | /[ e3w(I+1,J,K) + e3w(I,J,K) - e3w(I+1,J,K) ] 989 !! | | | 990 !! | | | 2. Horizontal interpolation: compute value at U/V point 991 !!K-1 | ptab(I,J,K-1) | | interpolation between zbis and ptab(I+1,J,K) 992 !! | . | | 993 !! | . | | Z= ( 0.5*zet2*zbis + 0.5*zet1*ptab(I+1,J,K) )/(0.5*zet2+0.5*zet1) 994 !! | . | | 995 !! ------------------------------------------ 996 !! | . | | 997 !! | . | | 998 !! | . | | 999 !!K | zbis...... U ..ptab(I+1,J,K) | 1000 !! | . | | 1001 !! | ptab(I,J,K) | | 1002 !! | |------------------| 1003 !! | | partials | 1004 !! | | steps | 1005 !! ------------------------------------------- 1006 !! <----zet1------><----zet2---------> 1007 !! 1008 !! 1009 !! ====> s-coordinate 1010 !! 1011 !! | | | 1. Compute distance between T1 and U points 1012 !! | | | Compute distance between T2 and U points 1013 !! | | ptab(I+1,J,K) | 1014 !! | | T2 | 2. Interpolation between T1 and T2 values at U point 1015 !! | | ^ | 1016 !! | | | zdep2 | 1017 !! | | | | 1018 !! | ^ U v | 1019 !! | | | | 1020 !! | | zdep1 | | 1021 !! | v | | 1022 !! | T1 | | 1023 !! | ptab(I,J,K) | | 1024 !! | | | 1025 !! | | | 1026 !! 1027 !! <----zet1--------><----zet2---------> 997 1028 !! 998 1029 !!---------------------------------------------------------------------- 999 1030 !*arguments 1000 INTEGER, INTENT(IN) :: ki, kj, kk 1001 CHARACTER(len=1), INTENT(IN) :: cd_point 1002 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(IN) :: ptab 1003 REAL(wp) :: interp 1031 INTEGER, INTENT(IN) :: ki, kj, kk ! coordinate of point 1032 CHARACTER(len=1), INTENT(IN) :: cd_point ! type of point (U, V) 1033 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(IN) :: ptab ! variable to compute at (ki, kj, kk ) 1034 REAL(wp) :: interp ! interpolated variable 1004 1035 1005 1036 !*local declations 1006 INTEGER :: ii1, ij1, ii2, ij2 1007 REAL(wp):: ze3w, zfse3, zmax1, zmax2, zbis 1037 INTEGER :: ii1, ij1, ii2, ij2 ! local integer 1038 REAL(wp):: ze3t, zfse3, zwgt1, zwgt2, zbis, zdepu ! local real 1039 REAL(wp):: zet1, zet2 ! weight for interpolation 1040 REAL(wp):: zdep1,zdep2 ! differences of depth 1008 1041 !!---------------------------------------------------------------------- 1009 1042 … … 1011 1044 ii1 = ki ; ij1 = kj 1012 1045 ii2 = ki+1 ; ij2 = kj 1046 1047 zet1=e1t(ii1,ij1) 1048 zet2=e1t(ii2,ij2) 1049 1013 1050 ELSE ! cd_point=='V' 1014 1051 ii1 = ki ; ij1 = kj 1015 1052 ii2 = ki ; ij2 = kj+1 1053 1054 zet1=e2t(ii1,ij1) 1055 zet2=e2t(ii2,ij2) 1056 1016 1057 ENDIF 1017 1058 1059 IF( ln_sco )THEN ! s-coordinate case 1060 1061 zdepu = ( fsdept(ii1,ij1,kk) + fsdept(ii2,ij2,kk) ) /2 1062 zdep1 = fsdept(ii1,ij1,kk) - zdepu 1063 zdep2 = fsdept(ii2,ij2,kk) - zdepu 1064 1065 zwgt1 = SQRT( ( 0.5 * zet1 ) * ( 0.5 * zet1 ) + ( zdep1 * zdep1 ) ) 1066 zwgt2 = SQRT( ( 0.5 * zet2 ) * ( 0.5 * zet2 ) + ( zdep2 * zdep2 ) ) 1067 1068 ! result 1069 interp = umask(ii1,ij1,kk) * ( zwgt2 * ptab(ii1,ij1,kk) + zwgt1 * ptab(ii1,ij1,kk) ) / ( zwgt2 + zwgt1 ) 1070 1071 1072 ELSE ! full step or partial step case 1073 1018 1074 #if defined key_vvl 1019 1075 1020 ze3w = fsve3w(ii2,ij2,kk) - fsve3w(ii1,ij1,kk)1021 1022 zfse3 = fsve3w(ii1,ij1,kk) * ( 1 + sshn(ii2,ij2) * mut(ii2,ij2,kk) )1023 zmax1 = ( fse3w(ii2,ij2,kk) - zfse3 )/ fse3w(ii2,ij2,kk)1076 ze3w = fsve3w(ii2,ij2,kk) - fsve3w(ii1,ij1,kk) 1077 1078 zfse3 = fsve3w(ii1,ij1,kk) * ( 1 + sshn(ii2,ij2) * mut(ii2,ij2,kk) ) 1079 zwgt1 = ( fse3w(ii2,ij2,kk) - zfse3 )/ fse3w(ii2,ij2,kk) 1024 1080 1025 zfse3 = fsve3w(ii2,ij2,kk) * ( 1 + sshn(ii1,ij1) * mut(ii1,ij1,kk) )1026 zmax2 = ( fse3w(ii1,ij1,kk) - zfse3 )/ fse3w(ii1,ij1,kk)1081 zfse3 = fsve3w(ii2,ij2,kk) * ( 1 + sshn(ii1,ij1) * mut(ii1,ij1,kk) ) 1082 zwgt2 = ( fse3w(ii1,ij1,kk) - zfse3 )/ fse3w(ii1,ij1,kk) 1027 1083 1028 1084 #else 1029 1085 1030 ze3w = fse3w(ii2,ij2,kk)-fse3w(ii1,ij1,kk)1031 zmax1 = ze3w/ fse3w(ii2,ij2,kk)1032 zmax2 = -ze3w/ fse3w(ii1,ij1,kk)1086 ze3t = fse3t(ii2,ij2,kk)-fse3t(ii1,ij1,kk) 1087 zwgt1 = ( fse3w(ii2,ij2,kk) - fse3w(ii1,ij1,kk) ) / fse3w(ii2,ij2,kk) 1088 zwgt2 = ( fse3w(ii1,ij1,kk) - fse3w(ii2,ij2,kk) ) / fse3w(ii1,ij1,kk) 1033 1089 1034 1090 #endif 1035 1091 1036 IF(kk .NE. 1)THEN 1037 1038 IF( ze3w >= 0. )THEN 1039 !zbis 1040 zbis = ptab(ii2,ij2,kk) + zmax1 * ( ptab(ii2,ij2,kk-1) - ptab(ii2,ij2,kk) ) 1041 ! result 1042 interp = umask(ii1,ij1,kk) * 0.5*( ptab(ii1,ij1,kk) + zbis ) 1092 IF(kk .NE. 1)THEN 1093 1094 IF( ze3t >= 0. )THEN 1095 !zbis 1096 zbis = ptab(ii2,ij2,kk) + zwgt1 * ( ptab(ii2,ij2,kk-1) - ptab(ii2,ij2,kk) ) 1097 ! result 1098 interp = umask(ii1,ij1,kk) * ( zet2 * ptab(ii1,ij1,kk) + zet1 * zbis )/( zet1 + zet2 ) 1099 ELSE 1100 !zbis 1101 zbis = ptab(ii1,ij1,kk) + zwgt2 * ( ptab(ii1,ij1,kk-1) - ptab(ii1,ij2,kk) ) 1102 ! result 1103 interp = umask(ii1,ij1,kk) * ( zet2 * zbis + zet1 * ptab(ii2,ij2,kk) )/( zet1 + zet2 ) 1104 ENDIF 1105 1043 1106 ELSE 1044 !zbis 1045 zbis = ptab(ii1,ij1,kk) + zmax2 * ( ptab(ii1,ij1,kk-1) - ptab(ii1,ij2,kk) ) 1046 ! result 1047 interp = umask(ii1,ij1,kk) * 0.5*( zbis + ptab(ii2,ij2,kk) ) 1048 ENDIF 1049 1050 ELSE 1051 interp = umask(ii1,ij1,kk) * 0.5*( ptab(ii1,ij1,kk) + ptab(ii2,ij2,kk) ) 1107 interp = umask(ii1,ij1,kk) * ( zet2 * ptab(ii1,ij1,kk) + zet1 * ptab(ii2,ij2,kk) )/( zet1 + zet2 ) 1108 ENDIF 1109 1052 1110 ENDIF 1111 1053 1112 1054 1113 END FUNCTION interp
Note: See TracChangeset
for help on using the changeset viewer.