Changeset 10237
- Timestamp:
- 2018-10-26T17:21:47+02:00 (6 years ago)
- Location:
- branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90
r8058 r10237 59 59 PRIVATE transport 60 60 PRIVATE dia_dct_wri 61 PRIVATE dia_dct_wri_NOOS 61 62 62 63 #include "domzgr_substitute.h90" … … 66 67 67 68 !! * Module variables 68 INTEGER :: nn_dct ! Frequency of computation 69 INTEGER :: nn_dctwri ! Frequency of output 70 INTEGER :: nn_secdebug ! Number of the section to debug 69 INTEGER :: nn_dct = 1 ! Frequency of computation 70 INTEGER :: nn_dctwri = 1 ! Frequency of output 71 INTEGER :: nn_secdebug = 0 ! Number of the section to debug 72 INTEGER :: nn_dct_h = 1 ! Frequency of computation for NOOS hourly files 73 INTEGER :: nn_dctwri_h = 1 ! Frequency of output for NOOS hourly files 71 74 72 INTEGER, PARAMETER :: nb_class_max = 1 073 INTEGER, PARAMETER :: nb_sec_max = 15074 INTEGER, PARAMETER :: nb_point_max = 200075 INTEGER, PARAMETER :: nb_type_class = 1076 INTEGER, PARAMETER :: nb_3d_vars = 377 INTEGER, PARAMETER :: nb_2d_vars = 2 75 INTEGER, PARAMETER :: nb_class_max = 12 ! maximum number of classes, i.e. depth levels or density classes 76 INTEGER, PARAMETER :: nb_sec_max = 30 ! maximum number of sections 77 INTEGER, PARAMETER :: nb_point_max = 300 ! maximum number of points in a single section 78 INTEGER, PARAMETER :: nb_type_class = 14 ! types of calculations, i.e. pos transport, neg transport, heat transport, salt transport 79 INTEGER, PARAMETER :: nb_3d_vars = 5 80 INTEGER, PARAMETER :: nb_2d_vars = 2 78 81 INTEGER :: nb_sec 79 82 … … 102 105 zlay ! level classes (99 if you don't want) 103 106 REAL(wp), DIMENSION(nb_type_class,nb_class_max) :: transport ! transport output 107 REAL(wp), DIMENSION(nb_type_class,nb_class_max) :: transport_h ! transport output 104 108 REAL(wp) :: slopeSection ! slope of the section 105 109 INTEGER :: nb_point ! number of points in the section … … 111 115 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: transports_3d 112 116 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: transports_2d 117 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: transports_3d_h 118 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: transports_2d_h 119 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: z_hr_output 113 120 114 121 !! $Id$ … … 120 127 !! *** FUNCTION diadct_alloc *** 121 128 !!---------------------------------------------------------------------- 122 INTEGER :: ierr( 2)129 INTEGER :: ierr(5) 123 130 !!---------------------------------------------------------------------- 124 131 125 ALLOCATE(transports_3d(nb_3d_vars,nb_sec_max,nb_point_max,jpk), STAT=ierr(1) ) 126 ALLOCATE(transports_2d(nb_2d_vars,nb_sec_max,nb_point_max) , STAT=ierr(2) ) 132 ALLOCATE(transports_3d(nb_3d_vars,nb_sec_max,nb_point_max,jpk) , STAT=ierr(1) ) 133 ALLOCATE(transports_2d(nb_2d_vars,nb_sec_max,nb_point_max) , STAT=ierr(2) ) 134 ALLOCATE(transports_3d_h(nb_3d_vars,nb_sec_max,nb_point_max,jpk), STAT=ierr(3) ) 135 ALLOCATE(transports_2d_h(nb_2d_vars,nb_sec_max,nb_point_max) , STAT=ierr(4) ) 136 ALLOCATE(z_hr_output(nb_sec_max,168,nb_class_max) , STAT=ierr(5) ) ! 168 = 24 * 7days 127 137 128 138 diadct_alloc = MAXVAL( ierr ) … … 153 163 IF(lwm) WRITE ( numond, namdct ) 154 164 165 IF( ln_NOOS ) THEN 166 nn_dct=3600./rdt ! hard coded for NOOS transects, to give 25 hour means 167 nn_dctwri=86400./rdt 168 nn_dct_h=1 ! hard coded for NOOS transects, to give hourly data 169 nn_dctwri_h=3600./rdt 170 ENDIF 171 155 172 IF( lwp ) THEN 156 173 WRITE(numout,*) " " 157 174 WRITE(numout,*) "diadct_init: compute transports through sections " 158 175 WRITE(numout,*) "~~~~~~~~~~~~~~~~~~~~~" 159 WRITE(numout,*) " Frequency of computation: nn_dct = ",nn_dct 160 WRITE(numout,*) " Frequency of write: nn_dctwri = ",nn_dctwri 176 IF( ln_NOOS ) THEN 177 WRITE(numout,*) " Frequency of computation hard coded to be every hour: nn_dct = ",nn_dct 178 WRITE(numout,*) " Frequency of write hard coded to average 25 instantaneous hour values: nn_dctwri = ",nn_dctwri 179 WRITE(numout,*) " Frequency of hourly computation hard coded to be every timestep: nn_dct_h = ",nn_dct_h 180 WRITE(numout,*) " Frequency of hourly write hard coded to every hour: nn_dctwri_h = ",nn_dctwri_h 181 ELSE 182 WRITE(numout,*) " Frequency of computation: nn_dct = ",nn_dct 183 WRITE(numout,*) " Frequency of write: nn_dctwri = ",nn_dctwri 184 ENDIF 161 185 162 186 IF ( nn_secdebug .GE. 1 .AND. nn_secdebug .LE. nb_sec_max )THEN … … 177 201 !open output file 178 202 IF( lwm ) THEN 179 CALL ctl_opn( numdct_vol, 'volume_transport', 'NEW', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 180 CALL ctl_opn( numdct_heat, 'heat_transport' , 'NEW', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 181 CALL ctl_opn( numdct_salt, 'salt_transport' , 'NEW', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 203 IF( ln_NOOS ) THEN 204 CALL ctl_opn( numdct_NOOS ,'NOOS_transport' , 'NEW', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 205 CALL ctl_opn( numdct_NOOS_h,'NOOS_transport_h', 'NEW', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 206 ELSE 207 CALL ctl_opn( numdct_vol, 'volume_transport', 'NEW', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 208 CALL ctl_opn( numdct_heat, 'heat_transport' , 'NEW', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 209 CALL ctl_opn( numdct_salt, 'salt_transport' , 'NEW', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 210 ENDIF 182 211 ENDIF 183 212 184 213 ! Initialise arrays to zero 185 transports_3d(:,:,:,:)=0.0 186 transports_2d(:,:,:) =0.0 214 transports_3d(:,:,:,:) =0.0 215 transports_2d(:,:,:) =0.0 216 transports_3d_h(:,:,:,:)=0._wp 217 transports_2d_h(:,:,:) =0._wp 218 z_hr_output(:,:,:) =0._wp 187 219 188 220 IF( nn_timing == 1 ) CALL timing_stop('dia_dct_init') … … 229 261 CALL wrk_alloc( nb_sec_max,nb_type_class,nb_class_max , zsum ) 230 262 ENDIF 231 263 lldebug=.TRUE. 232 264 ! Initialise arrays 233 265 zwork(:) = 0.0 … … 243 275 244 276 ! Compute transport and write only at nn_dctwri 245 IF( MOD(kt,nn_dct)==0 ) THEN 277 IF( MOD(kt,nn_dct)==0 .or. & ! compute transport every nn_dct time steps 278 (ln_NOOS .and. kt==nn_it000 ) ) THEN ! also include first time step when calculating NOOS 25 hour averages 246 279 247 280 DO jsec=1,nb_sec 248 281 249 282 !debug this section computing ? 250 lldebug=.FALSE.251 283 IF( (jsec==nn_secdebug .OR. nn_secdebug==-1) .AND. kt==nit000+nn_dct-1 .AND. lwp ) lldebug=.TRUE. 252 284 … … 271 303 !Sum on all procs 272 304 IF( lk_mpp )THEN 305 zsum(:,:,:)=0.0_wp 273 306 ish(1) = nb_sec_max*nb_type_class*nb_class_max 274 307 ish2 = (/nb_sec_max,nb_type_class,nb_class_max/) … … 283 316 DO jsec=1,nb_sec 284 317 285 IF( lwm )CALL dia_dct_wri(kt,jsec,secs(jsec)) 318 IF( lwm .and. .not. ln_NOOS )CALL dia_dct_wri(kt,jsec,secs(jsec)) 319 IF( lwm .and. ln_NOOS )CALL dia_dct_wri_NOOS(kt,jsec,secs(jsec)) ! use NOOS specific formatting 286 320 287 321 !nullify transports values after writing 288 322 transports_3d(:,jsec,:,:)=0. 289 323 transports_2d(:,jsec,: )=0. 290 secs(jsec)%transport(:,:)=0. 324 secs(jsec)%transport(:,:)=0. 325 IF ( ln_NOOS ) CALL transport(secs(jsec),lldebug,jsec) ! reinitialise for next 25 hour instantaneous average (overlapping values) 291 326 292 327 ENDDO … … 295 330 296 331 ENDIF 332 333 IF ( MOD(kt,nn_dct_h)==0 ) THEN ! compute transport every nn_dct_h time steps 334 335 DO jsec=1,nb_sec 336 337 !lldebug=.FALSE. 338 IF( (jsec==nn_secdebug .OR. nn_secdebug==-1) .AND. kt==nit000+nn_dct_h-1 .AND. lwp ) lldebug=.TRUE. 339 340 !Compute transport through section 341 CALL transport_h(secs(jsec),lldebug,jsec) 342 343 ENDDO 344 345 IF( MOD(kt,nn_dctwri_h)==0 )THEN 346 347 IF( lwp .AND. kt==nit000+nn_dctwri_h-1 )WRITE(numout,*)" diadct: average and write hourly files at kt = ",kt 348 349 !! divide arrays by nn_dctwri/nn_dct to obtain average 350 transports_3d_h(:,:,:,:)=transports_3d_h(:,:,:,:)/(nn_dctwri_h/nn_dct_h) 351 transports_2d_h(:,:,:) =transports_2d_h(:,:,:) /(nn_dctwri_h/nn_dct_h) 352 353 ! Sum over each class 354 DO jsec=1,nb_sec 355 CALL dia_dct_sum_h(secs(jsec),jsec) 356 ENDDO 357 358 !Sum on all procs 359 IF( lk_mpp )THEN 360 ish(1) = nb_sec_max*nb_type_class*nb_class_max 361 ish2 = (/nb_sec_max,nb_type_class,nb_class_max/) 362 DO jsec=1,nb_sec ; zsum(jsec,:,:) = secs(jsec)%transport_h(:,:) ; ENDDO 363 zwork(:)= RESHAPE(zsum(:,:,:), ish ) 364 CALL mpp_sum(zwork, ish(1)) 365 zsum(:,:,:)= RESHAPE(zwork,ish2) 366 DO jsec=1,nb_sec ; secs(jsec)%transport_h(:,:) = zsum(jsec,:,:) ; ENDDO 367 ENDIF 368 369 !Write the transport 370 DO jsec=1,nb_sec 371 372 IF( lwp .and. ln_NOOS )CALL dia_dct_wri_NOOS_h(kt/nn_dctwri_h,jsec,secs(jsec)) ! use NOOS specific formatting 373 374 !nullify transports values after writing 375 transports_3d_h(:,jsec,:,:)=0.0 376 transports_2d_h(:,jsec,:)=0.0 377 secs(jsec)%transport_h(:,:)=0. 378 IF ( ln_NOOS ) CALL transport_h(secs(jsec),lldebug,jsec) ! reinitialise for next 25 hour instantaneous average (overlapping values) 379 380 ENDDO 381 382 ENDIF 383 384 ENDIF 297 385 298 386 IF( lk_mpp )THEN … … 356 444 secs(jsec)%zlay = 99._wp 357 445 secs(jsec)%transport = 0._wp ; secs(jsec)%nb_point = 0 446 secs(jsec)%transport_h = 0._wp ; secs(jsec)%nb_point = 0 358 447 359 448 !read section's number / name / computing choices / classes / slopeSection / points number 360 449 !----------------------------------------------------------------------------------------- 361 450 READ(numdct_in,iostat=iost)isec 362 IF (iost .NE. 0 )EXIT !end of file 451 IF (iost .NE. 0 )EXIT !end of file 363 452 WRITE(cltmp,'(a,i4.4,a,i4.4)')'diadct: read sections : Problem of section number: isec= ',isec,' and jsec= ',jsec 364 IF( jsec .NE. isec ) 453 IF( jsec .NE. isec )CALL ctl_stop( cltmp ) 365 454 366 455 IF( lwp .AND. ( jsec==nn_secdebug .OR. nn_secdebug==-1 ) )WRITE(numout,*)"isec ",isec … … 380 469 READ(numdct_in)iptglo 381 470 471 IF ( ln_NOOS ) THEN 472 WRITE(numout,*) 'Section name = ',TRIM(secs(jsec)%name) 473 WRITE(numout,*) 'coordSec = ',secs(jsec)%coordSec 474 WRITE(numout,*) 'iptglo = ',iptglo 475 ENDIF 476 382 477 !debug 383 478 !----- … … 405 500 !read points'coordinates and directions 406 501 !-------------------------------------- 502 IF ( ln_NOOS ) THEN 503 WRITE(numout,*) 'Coords and direction read' 504 ENDIF 505 407 506 coordtemp(:) = POINT_SECTION(0,0) !list of points read 408 507 directemp(:) = 0 !value of directions of each points … … 422 521 ENDDO 423 522 ENDIF 424 523 425 524 !Now each proc selects only points that are in its domain: 426 525 !-------------------------------------------------------- … … 624 723 !!-------------------------------------------------------- 625 724 626 IF( ld_debug )WRITE(numout,*)' Compute transport'725 !!NIALL IF( ld_debug )WRITE(numout,*)' Compute transport' 627 726 628 727 !---------------------------! … … 710 809 SELECT CASE( sec%direction(jseg) ) 711 810 CASE(0,1) 712 ztn = interp(k%I,k%J,jk,'V', tsn(:,:,:,jp_tem))713 zsn = interp(k%I,k%J,jk,'V', tsn(:,:,:,jp_sal))714 zrhop = interp(k%I,k%J,jk,'V', rhop)715 zrhoi = interp(k%I,k%J,jk,'V', rhd*rau0+rau0)811 ztn = interp(k%I,k%J,jk,'V',0 ) 812 zsn = interp(k%I,k%J,jk,'V',1 ) 813 zrhop = interp(k%I,k%J,jk,'V',2 ) 814 zrhoi = interp(k%I,k%J,jk,'V',3 ) 716 815 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) * vmask(k%I,k%J,1) 717 816 CASE(2,3) 718 ztn = interp(k%I,k%J,jk,'U', tsn(:,:,:,jp_tem))719 zsn = interp(k%I,k%J,jk,'U', tsn(:,:,:,jp_sal))720 zrhop = interp(k%I,k%J,jk,'U', rhop)721 zrhoi = interp(k%I,k%J,jk,'U', rhd*rau0+rau0)817 ztn = interp(k%I,k%J,jk,'U',0 ) 818 zsn = interp(k%I,k%J,jk,'U',1 ) 819 zrhop = interp(k%I,k%J,jk,'U',2 ) 820 zrhoi = interp(k%I,k%J,jk,'U',3 ) 722 821 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) * umask(k%I,k%J,1) 723 822 END SELECT … … 754 853 transports_3d(2,jsec,jseg,jk) = transports_3d(2,jsec,jseg,jk) + zTnorm * ztn * zrhop * rcp 755 854 transports_3d(3,jsec,jseg,jk) = transports_3d(3,jsec,jseg,jk) + zTnorm * zsn * zrhop * 0.001 855 transports_3d(4,jsec,jseg,jk) = transports_3d(4,jsec,jseg,jk) + zTnorm * 4.e+3_wp * (ztn+273.15) * 1026._wp 856 transports_3d(5,jsec,jseg,jk) = transports_3d(5,jsec,jseg,jk) + zTnorm * 0.001 * zsn * 1026._wp 756 857 ENDIF 757 858 … … 809 910 END SUBROUTINE transport 810 911 912 SUBROUTINE transport_h(sec,ld_debug,jsec) 913 !!------------------------------------------------------------------------------------------- 914 !! *** ROUTINE hourly_transport *** 915 !! 916 !! Exactly as routine transport but for data summed at 917 !! each time step and averaged each hour 918 !! 919 !! Purpose :: Compute the transport for each point in a section 920 !! 921 !! Method :: Loop over each segment, and each vertical level and add the transport 922 !! Be aware : 923 !! One section is a sum of segments 924 !! One segment is defined by 2 consecutive points in sec%listPoint 925 !! All points of sec%listPoint are positioned on the F-point of the cell 926 !! 927 !! There are two loops: 928 !! loop on the segment between 2 nodes 929 !! loop on the level jk 930 !! 931 !! ** Output: Arrays containing the volume,density,salinity,temperature etc 932 !! transports for each point in a section, summed over each nn_dct. 933 !! 934 !!------------------------------------------------------------------------------------------- 935 !! * Arguments 936 TYPE(SECTION),INTENT(INOUT) :: sec 937 LOGICAL ,INTENT(IN) :: ld_debug 938 INTEGER ,INTENT(IN) :: jsec ! numeric identifier of section 939 940 !! * Local variables 941 INTEGER :: jk,jseg,jclass, &!loop on level/segment/classes 942 isgnu , isgnv ! 943 REAL(wp):: zumid , zvmid , &!U/V velocity on a cell segment 944 zumid_ice , zvmid_ice , &!U/V ice velocity 945 zTnorm !transport of velocity through one cell's sides 946 REAL(wp):: ztn, zsn, zrhoi, zrhop, zsshn, zfsdep ! temperature/salinity/ssh/potential density /depth at u/v point 947 948 TYPE(POINT_SECTION) :: k 949 !!-------------------------------------------------------- 950 951 !!NIALL IF( ld_debug )WRITE(numout,*)' Compute transport' 952 953 !---------------------------! 954 ! COMPUTE TRANSPORT ! 955 !---------------------------! 956 IF(sec%nb_point .NE. 0)THEN 957 958 !---------------------------------------------------------------------------------------------------- 959 !Compute sign for velocities: 960 ! 961 !convention: 962 ! non horizontal section: direction + is toward left hand of section 963 ! horizontal section: direction + is toward north of section 964 ! 965 ! 966 ! slopeSection < 0 slopeSection > 0 slopeSection=inf slopeSection=0 967 ! ---------------- ----------------- --------------- -------------- 968 ! 969 ! isgnv=1 direction + 970 ! ______ _____ ______ 971 ! | //| | | direction + 972 ! | isgnu=1 // | |isgnu=1 |isgnu=1 /|\ 973 ! |_______ // ______| \\ | ---\ | 974 ! | | isgnv=-1 \\ | | ---/ direction + ____________ 975 ! | | __\\| | 976 ! | | direction + | isgnv=1 977 ! 978 !---------------------------------------------------------------------------------------------------- 979 isgnu = 1 980 IF( sec%slopeSection .GT. 0 ) THEN ; isgnv = -1 981 ELSE ; isgnv = 1 982 ENDIF 983 984 IF( ld_debug )write(numout,*)"isgnu isgnv ",isgnu,isgnv 985 986 !--------------------------------------! 987 ! LOOP ON THE SEGMENT BETWEEN 2 NODES ! 988 !--------------------------------------! 989 DO jseg=1,MAX(sec%nb_point-1,0) 990 991 !------------------------------------------------------------------------------------------- 992 ! Select the appropriate coordinate for computing the velocity of the segment 993 ! 994 ! CASE(0) Case (2) 995 ! ------- -------- 996 ! listPoint(jseg) listPoint(jseg+1) listPoint(jseg) F(i,j) 997 ! F(i,j)----------V(i+1,j)-------F(i+1,j) | 998 ! | 999 ! | 1000 ! | 1001 ! Case (3) U(i,j) 1002 ! -------- | 1003 ! | 1004 ! listPoint(jseg+1) F(i,j+1) | 1005 ! | | 1006 ! | | 1007 ! | listPoint(jseg+1) F(i,j-1) 1008 ! | 1009 ! | 1010 ! U(i,j+1) 1011 ! | Case(1) 1012 ! | ------ 1013 ! | 1014 ! | listPoint(jseg+1) listPoint(jseg) 1015 ! | F(i-1,j)-----------V(i,j) -------f(jseg) 1016 ! listPoint(jseg) F(i,j) 1017 ! 1018 !------------------------------------------------------------------------------------------- 1019 1020 SELECT CASE( sec%direction(jseg) ) 1021 CASE(0) ; k = sec%listPoint(jseg) 1022 CASE(1) ; k = POINT_SECTION(sec%listPoint(jseg)%I+1,sec%listPoint(jseg)%J) 1023 CASE(2) ; k = sec%listPoint(jseg) 1024 CASE(3) ; k = POINT_SECTION(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J+1) 1025 END SELECT 1026 1027 !---------------------------| 1028 ! LOOP ON THE LEVEL | 1029 !---------------------------| 1030 !Sum of the transport on the vertical 1031 DO jk=1,mbathy(k%I,k%J) 1032 1033 ! compute temparature, salinity, insitu & potential density, ssh and depth at U/V point 1034 SELECT CASE( sec%direction(jseg) ) 1035 CASE(0,1) 1036 ztn = interp(k%I,k%J,jk,'V',0 ) 1037 zsn = interp(k%I,k%J,jk,'V',1) 1038 zrhop = interp(k%I,k%J,jk,'V',2) 1039 zrhoi = interp(k%I,k%J,jk,'V',3) 1040 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) * vmask(k%I,k%J,1) 1041 CASE(2,3) 1042 ztn = interp(k%I,k%J,jk,'U',0) 1043 zsn = interp(k%I,k%J,jk,'U',1) 1044 zrhop = interp(k%I,k%J,jk,'U',2) 1045 zrhoi = interp(k%I,k%J,jk,'U',3) 1046 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) * umask(k%I,k%J,1) 1047 END SELECT 1048 1049 zfsdep= fsdept(k%I,k%J,jk) 1050 1051 !compute velocity with the correct direction 1052 SELECT CASE( sec%direction(jseg) ) 1053 CASE(0,1) 1054 zumid=0. 1055 zvmid=isgnv*vn(k%I,k%J,jk)*vmask(k%I,k%J,jk) 1056 CASE(2,3) 1057 zumid=isgnu*un(k%I,k%J,jk)*umask(k%I,k%J,jk) 1058 zvmid=0. 1059 END SELECT 1060 1061 !zTnorm=transport through one cell; 1062 !velocity* cell's length * cell's thickness 1063 zTnorm=zumid*e2u(k%I,k%J)* fse3u(k%I,k%J,jk)+ & 1064 zvmid*e1v(k%I,k%J)* fse3v(k%I,k%J,jk) 1065 1066 #if ! defined key_vvl 1067 !add transport due to free surface 1068 IF( jk==1 )THEN 1069 zTnorm = zTnorm + zumid* e2u(k%I,k%J) * zsshn * umask(k%I,k%J,jk) + & 1070 zvmid* e1v(k%I,k%J) * zsshn * vmask(k%I,k%J,jk) 1071 ENDIF 1072 #endif 1073 !COMPUTE TRANSPORT 1074 1075 transports_3d_h(1,jsec,jseg,jk) = transports_3d_h(1,jsec,jseg,jk) + zTnorm 1076 1077 IF ( sec%llstrpond ) THEN 1078 transports_3d_h(2,jsec,jseg,jk) = transports_3d_h(2,jsec,jseg,jk) + zTnorm * zrhoi 1079 transports_3d_h(3,jsec,jseg,jk) = transports_3d_h(3,jsec,jseg,jk) + zTnorm * zrhop 1080 transports_3d_h(4,jsec,jseg,jk) = transports_3d_h(4,jsec,jseg,jk) + zTnorm * 4.e+3_wp * (ztn+273.15) * 1026._wp 1081 transports_3d_h(5,jsec,jseg,jk) = transports_3d_h(5,jsec,jseg,jk) + zTnorm * 0.001 * zsn * 1026._wp 1082 ENDIF 1083 1084 ENDDO !end of loop on the level 1085 1086 #if defined key_lim2 || defined key_lim3 1087 1088 !ICE CASE 1089 !------------ 1090 IF( sec%ll_ice_section )THEN 1091 SELECT CASE (sec%direction(jseg)) 1092 CASE(0) 1093 zumid_ice = 0 1094 zvmid_ice = isgnv*0.5*(v_ice(k%I,k%J+1)+v_ice(k%I+1,k%J+1)) 1095 CASE(1) 1096 zumid_ice = 0 1097 zvmid_ice = isgnv*0.5*(v_ice(k%I,k%J+1)+v_ice(k%I+1,k%J+1)) 1098 CASE(2) 1099 zvmid_ice = 0 1100 zumid_ice = isgnu*0.5*(u_ice(k%I+1,k%J)+u_ice(k%I+1,k%J+1)) 1101 CASE(3) 1102 zvmid_ice = 0 1103 zumid_ice = isgnu*0.5*(u_ice(k%I+1,k%J)+u_ice(k%I+1,k%J+1)) 1104 END SELECT 1105 1106 zTnorm=zumid_ice*e2u(k%I,k%J)+zvmid_ice*e1v(k%I,k%J) 1107 1108 transports_2d_h(1,jsec,jseg) = transports_2d_h(1,jsec,jseg) + (zTnorm)* & 1109 (1.0 - frld(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) & 1110 *(hsnif(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J) + & 1111 hicif(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) & 1112 +zice_vol_pos 1113 transports_2d_h(2,jsec,jseg) = transports_2d_h(2,jsec,jseg) + (zTnorm)* & 1114 (1.0 - frld(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J)) & 1115 +zice_surf_pos 1116 1117 ENDIF !end of ice case 1118 #endif 1119 1120 ENDDO !end of loop on the segment 1121 1122 ENDIF !end of sec%nb_point =0 case 1123 ! 1124 END SUBROUTINE transport_h 1125 811 1126 SUBROUTINE dia_dct_sum(sec,jsec) 812 1127 !!------------------------------------------------------------- … … 890 1205 SELECT CASE( sec%direction(jseg) ) 891 1206 CASE(0,1) 892 ztn = interp(k%I,k%J,jk,'V', tsn(:,:,:,jp_tem))893 zsn = interp(k%I,k%J,jk,'V', tsn(:,:,:,jp_sal))894 zrhop = interp(k%I,k%J,jk,'V', rhop)895 zrhoi = interp(k%I,k%J,jk,'V', rhd*rau0+rau0)1207 ztn = interp(k%I,k%J,jk,'V',0 ) 1208 zsn = interp(k%I,k%J,jk,'V',1 ) 1209 zrhop = interp(k%I,k%J,jk,'V',2) 1210 zrhoi = interp(k%I,k%J,jk,'V',3) 896 1211 897 1212 CASE(2,3) 898 ztn = interp(k%I,k%J,jk,'U', tsn(:,:,:,jp_tem))899 zsn = interp(k%I,k%J,jk,'U', tsn(:,:,:,jp_sal))900 zrhop = interp(k%I,k%J,jk,'U', rhop)901 zrhoi = interp(k%I,k%J,jk,'U', rhd*rau0+rau0)1213 ztn = interp(k%I,k%J,jk,'U',0 ) 1214 zsn = interp(k%I,k%J,jk,'U',1 ) 1215 zrhop = interp(k%I,k%J,jk,'U',2 ) 1216 zrhoi = interp(k%I,k%J,jk,'U',3 ) 902 1217 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) * umask(k%I,k%J,1) 903 1218 END SELECT … … 957 1272 sec%transport(6,jclass) = sec%transport(6,jclass)+transports_3d(3,jsec,jseg,jk) 958 1273 ENDIF 1274 1275 IF ( transports_3d(4,jsec,jseg,jk) .GE. 0.0 ) THEN 1276 sec%transport(7,jclass) = sec%transport(7,jclass)+transports_3d(4,jsec,jseg,jk) 1277 ELSE 1278 sec%transport(8,jclass) = sec%transport(8,jclass)+transports_3d(4,jsec,jseg,jk) 1279 ENDIF 1280 1281 IF ( transports_3d(5,jsec,jseg,jk) .GE. 0.0 ) THEN 1282 sec%transport( 9,jclass) = sec%transport( 9,jclass)+transports_3d(5,jsec,jseg,jk) 1283 ELSE 1284 sec%transport(10,jclass) = sec%transport(10,jclass)+transports_3d(5,jsec,jseg,jk) 1285 ENDIF 959 1286 960 1287 ELSE … … 963 1290 sec%transport( 5,jclass) = 0._wp 964 1291 sec%transport( 6,jclass) = 0._wp 1292 sec%transport( 7,jclass) = 0._wp 1293 sec%transport( 8,jclass) = 0._wp 1294 sec%transport( 9,jclass) = 0._wp 1295 sec%transport(10,jclass) = 0._wp 965 1296 ENDIF 966 1297 … … 977 1308 978 1309 IF ( transports_2d(1,jsec,jseg) .GE. 0.0 ) THEN 979 sec%transport( 7,1) = sec%transport( 7,1)+transports_2d(1,jsec,jseg)*1.E-61310 sec%transport(11,1) = sec%transport(11,1)+transports_2d(1,jsec,jseg)*1.E-6 980 1311 ELSE 981 sec%transport( 8,1) = sec%transport( 8,1)+transports_2d(1,jsec,jseg)*1.E-61312 sec%transport(12,1) = sec%transport(12,1)+transports_2d(1,jsec,jseg)*1.E-6 982 1313 ENDIF 983 1314 984 1315 IF ( transports_2d(3,jsec,jseg) .GE. 0.0 ) THEN 985 sec%transport( 9,1) = sec%transport( 9,1)+transports_2d(2,jsec,jseg)*1.E-61316 sec%transport(13,1) = sec%transport(13,1)+transports_2d(2,jsec,jseg)*1.E-6 986 1317 ELSE 987 sec%transport(1 0,1) = sec%transport(10,1)+transports_2d(2,jsec,jseg)*1.E-61318 sec%transport(14,1) = sec%transport(14,1)+transports_2d(2,jsec,jseg)*1.E-6 988 1319 ENDIF 989 1320 … … 995 1326 ELSE !if sec%nb_point =0 996 1327 sec%transport(1:2,:)=0. 997 IF (sec%llstrpond) sec%transport(3: 6,:)=0.998 IF (sec%ll_ice_section) sec%transport( 7:10,:)=0.1328 IF (sec%llstrpond) sec%transport(3:10,:)=0. 1329 IF (sec%ll_ice_section) sec%transport(11:14,:)=0. 999 1330 ENDIF !end of sec%nb_point =0 case 1000 1331 1001 1332 END SUBROUTINE dia_dct_sum 1333 1334 SUBROUTINE dia_dct_sum_h(sec,jsec) 1335 !!------------------------------------------------------------- 1336 !! Exactly as dia_dct_sum but for hourly files containing data summed at each time step 1337 !! 1338 !! Purpose: Average the transport over nn_dctwri time steps 1339 !! and sum over the density/salinity/temperature/depth classes 1340 !! 1341 !! Method: 1342 !! Sum over relevant grid cells to obtain values 1343 !! for each 1344 !! There are several loops: 1345 !! loop on the segment between 2 nodes 1346 !! loop on the level jk 1347 !! loop on the density/temperature/salinity/level classes 1348 !! test on the density/temperature/salinity/level 1349 !! 1350 !! ** Method :Transport through a given section is equal to the sum of transports 1351 !! computed on each proc. 1352 !! On each proc,transport is equal to the sum of transport computed through 1353 !! segments linking each point of sec%listPoint with the next one. 1354 !! 1355 !!------------------------------------------------------------- 1356 !! * arguments 1357 TYPE(SECTION),INTENT(INOUT) :: sec 1358 INTEGER ,INTENT(IN) :: jsec ! numeric identifier of section 1359 1360 TYPE(POINT_SECTION) :: k 1361 INTEGER :: jk,jseg,jclass !loop on level/segment/classes 1362 REAL(wp) :: ztn, zsn, zrhoi, zrhop, zsshn, zfsdep ! temperature/salinity/ssh/potential density /depth at u/v point 1363 !!------------------------------------------------------------- 1364 1365 !! Sum the relevant segments to obtain values for each class 1366 IF(sec%nb_point .NE. 0)THEN 1367 1368 !--------------------------------------! 1369 ! LOOP ON THE SEGMENT BETWEEN 2 NODES ! 1370 !--------------------------------------! 1371 DO jseg=1,MAX(sec%nb_point-1,0) 1372 1373 !------------------------------------------------------------------------------------------- 1374 ! Select the appropriate coordinate for computing the velocity of the segment 1375 ! 1376 ! CASE(0) Case (2) 1377 ! ------- -------- 1378 ! listPoint(jseg) listPoint(jseg+1) listPoint(jseg) F(i,j) 1379 ! F(i,j)----------V(i+1,j)-------F(i+1,j) | 1380 ! | 1381 ! | 1382 ! | 1383 ! Case (3) U(i,j) 1384 ! -------- | 1385 ! | 1386 ! listPoint(jseg+1) F(i,j+1) | 1387 ! | | 1388 ! | | 1389 ! | listPoint(jseg+1) F(i,j-1) 1390 ! | 1391 ! | 1392 ! U(i,j+1) 1393 ! | Case(1) 1394 ! | ------ 1395 ! | 1396 ! | listPoint(jseg+1) listPoint(jseg) 1397 ! | F(i-1,j)-----------V(i,j) -------f(jseg) 1398 ! listPoint(jseg) F(i,j) 1399 ! 1400 !------------------------------------------------------------------------------------------- 1401 1402 SELECT CASE( sec%direction(jseg) ) 1403 CASE(0) ; k = sec%listPoint(jseg) 1404 CASE(1) ; k = POINT_SECTION(sec%listPoint(jseg)%I+1,sec%listPoint(jseg)%J) 1405 CASE(2) ; k = sec%listPoint(jseg) 1406 CASE(3) ; k = POINT_SECTION(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J+1) 1407 END SELECT 1408 1409 !---------------------------| 1410 ! LOOP ON THE LEVEL | 1411 !---------------------------| 1412 !Sum of the transport on the vertical 1413 DO jk=1,mbathy(k%I,k%J) 1414 1415 ! compute temperature, salinity, insitu & potential density, ssh and depth at U/V point 1416 SELECT CASE( sec%direction(jseg) ) 1417 CASE(0,1) 1418 ztn = interp(k%I,k%J,jk,'V',0 ) 1419 zsn = interp(k%I,k%J,jk,'V',1 ) 1420 zrhop = interp(k%I,k%J,jk,'V',2 ) 1421 zrhoi = interp(k%I,k%J,jk,'V',3 ) 1422 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) * vmask(k%I,k%J,1) 1423 CASE(2,3) 1424 ztn = interp(k%I,k%J,jk,'U',0 ) 1425 zsn = interp(k%I,k%J,jk,'U',1 ) 1426 zrhop = interp(k%I,k%J,jk,'U',2 ) 1427 zrhoi = interp(k%I,k%J,jk,'U',3 ) 1428 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) * umask(k%I,k%J,1) 1429 END SELECT 1430 1431 zfsdep= fsdept(k%I,k%J,jk) 1432 1433 !------------------------------- 1434 ! LOOP ON THE DENSITY CLASSES | 1435 !------------------------------- 1436 !The computation is made for each density/heat/salt/... class 1437 DO jclass=1,MAX(1,sec%nb_class-1) 1438 1439 !----------------------------------------------! 1440 !TEST ON THE DENSITY/SALINITY/TEMPERATURE/LEVEL! 1441 !----------------------------------------------! 1442 1443 IF ( ( & 1444 ((( zrhop .GE. (sec%zsigp(jclass)+1000. )) .AND. & 1445 ( zrhop .LE. (sec%zsigp(jclass+1)+1000. ))) .OR. & 1446 ( sec%zsigp(jclass) .EQ. 99.)) .AND. & 1447 1448 ((( zrhoi .GE. (sec%zsigi(jclass) + 1000. )) .AND. & 1449 ( zrhoi .LE. (sec%zsigi(jclass+1)+1000. ))) .OR. & 1450 ( sec%zsigi(jclass) .EQ. 99.)) .AND. & 1451 1452 ((( zsn .GT. sec%zsal(jclass)) .AND. & 1453 ( zsn .LE. sec%zsal(jclass+1))) .OR. & 1454 ( sec%zsal(jclass) .EQ. 99.)) .AND. & 1455 1456 ((( ztn .GE. sec%ztem(jclass)) .AND. & 1457 ( ztn .LE. sec%ztem(jclass+1))) .OR. & 1458 ( sec%ztem(jclass) .EQ.99.)) .AND. & 1459 1460 ((( zfsdep .GE. sec%zlay(jclass)) .AND. & 1461 ( zfsdep .LE. sec%zlay(jclass+1))) .OR. & 1462 ( sec%zlay(jclass) .EQ. 99. )) & 1463 )) THEN 1464 1465 !SUM THE TRANSPORTS FOR EACH CLASSES FOR THE POSITIVE AND NEGATIVE DIRECTIONS 1466 !---------------------------------------------------------------------------- 1467 IF (transports_3d_h(1,jsec,jseg,jk) .GE. 0.0) THEN 1468 sec%transport_h(1,jclass) = sec%transport_h(1,jclass)+transports_3d_h(1,jsec,jseg,jk) 1469 ELSE 1470 sec%transport_h(2,jclass) = sec%transport_h(2,jclass)+transports_3d_h(1,jsec,jseg,jk) 1471 ENDIF 1472 IF( sec%llstrpond )THEN 1473 1474 IF ( transports_3d_h(2,jsec,jseg,jk) .GE. 0.0 ) THEN 1475 sec%transport_h(3,jclass) = sec%transport_h(3,jclass)+transports_3d_h(2,jsec,jseg,jk) 1476 ELSE 1477 sec%transport_h(4,jclass) = sec%transport_h(4,jclass)+transports_3d_h(2,jsec,jseg,jk) 1478 ENDIF 1479 1480 IF ( transports_3d_h(3,jsec,jseg,jk) .GE. 0.0 ) THEN 1481 sec%transport_h(5,jclass) = sec%transport_h(5,jclass)+transports_3d_h(3,jsec,jseg,jk) 1482 ELSE 1483 sec%transport_h(6,jclass) = sec%transport_h(6,jclass)+transports_3d_h(3,jsec,jseg,jk) 1484 ENDIF 1485 1486 IF ( transports_3d_h(4,jsec,jseg,jk) .GE. 0.0 ) THEN 1487 sec%transport_h(7,jclass) = sec%transport_h(7,jclass)+transports_3d_h(4,jsec,jseg,jk) 1488 ELSE 1489 sec%transport_h(8,jclass) = sec%transport_h(8,jclass)+transports_3d_h(4,jsec,jseg,jk) 1490 ENDIF 1491 1492 IF ( transports_3d_h(5,jsec,jseg,jk) .GE. 0.0 ) THEN 1493 sec%transport_h( 9,jclass) = sec%transport_h( 9,jclass)+transports_3d_h(5,jsec,jseg,jk) 1494 ELSE 1495 sec%transport_h(10,jclass) = sec%transport_h(10,jclass)+transports_3d_h(5,jsec,jseg,jk) 1496 ENDIF 1497 1498 ELSE 1499 sec%transport_h( 3,jclass) = 0._wp 1500 sec%transport_h( 4,jclass) = 0._wp 1501 sec%transport_h( 5,jclass) = 0._wp 1502 sec%transport_h( 6,jclass) = 0._wp 1503 sec%transport_h( 7,jclass) = 0._wp 1504 sec%transport_h( 8,jclass) = 0._wp 1505 sec%transport_h( 9,jclass) = 0._wp 1506 sec%transport_h(10,jclass) = 0._wp 1507 ENDIF 1508 1509 ENDIF ! end of test if point is in class 1510 1511 ENDDO ! end of loop on the classes 1512 1513 ENDDO ! loop over jk 1514 1515 #if defined key_lim2 || defined key_lim3 1516 1517 !ICE CASE 1518 IF( sec%ll_ice_section )THEN 1519 1520 IF ( transports_2d_h(1,jsec,jseg) .GE. 0.0 ) THEN 1521 sec%transport_h(11,1) = sec%transport_h(11,1)+transports_2d_h(1,jsec,jseg) 1522 ELSE 1523 sec%transport_h(12,1) = sec%transport_h(12,1)+transports_2d_h(1,jsec,jseg) 1524 ENDIF 1525 1526 IF ( transports_2d_h(3,jsec,jseg) .GE. 0.0 ) THEN 1527 sec%transport_h(13,1) = sec%transport_h(13,1)+transports_2d_h(2,jsec,jseg) 1528 ELSE 1529 sec%transport_h(14,1) = sec%transport_h(14,1)+transports_2d_h(2,jsec,jseg) 1530 ENDIF 1531 1532 ENDIF !end of ice case 1533 #endif 1534 1535 ENDDO !end of loop on the segment 1536 1537 ELSE !if sec%nb_point =0 1538 sec%transport_h(1:2,:)=0. 1539 IF (sec%llstrpond) sec%transport_h(3:10,:)=0. 1540 IF (sec%ll_ice_section) sec%transport_h( 11:14,:)=0. 1541 ENDIF !end of sec%nb_point =0 case 1542 1543 END SUBROUTINE dia_dct_sum_h 1002 1544 1545 SUBROUTINE dia_dct_wri_NOOS(kt,ksec,sec) 1546 !!------------------------------------------------------------- 1547 !! Write transport output in numdct using NOOS formatting 1548 !! 1549 !! Purpose: Write transports in ascii files 1550 !! 1551 !! Method: 1552 !! 1. Write volume transports in "volume_transport" 1553 !! Unit: Sv : area * Velocity / 1.e6 1554 !! 1555 !! 2. Write heat transports in "heat_transport" 1556 !! Unit: Peta W : area * Velocity * T * rhau * Cp / 1.e15 1557 !! 1558 !! 3. Write salt transports in "salt_transport" 1559 !! Unit: 10^9 g m^3 / s : area * Velocity * S / 1.e6 1560 !! 1561 !!------------------------------------------------------------- 1562 !!arguments 1563 INTEGER, INTENT(IN) :: kt ! time-step 1564 TYPE(SECTION), INTENT(INOUT) :: sec ! section to write 1565 INTEGER ,INTENT(IN) :: ksec ! section number 1566 1567 !!local declarations 1568 INTEGER :: jclass,ji ! Dummy loop 1569 CHARACTER(len=2) :: classe ! Classname 1570 REAL(wp) :: zbnd1,zbnd2 ! Class bounds 1571 REAL(wp) :: zslope ! section's slope coeff 1572 ! 1573 REAL(wp), POINTER, DIMENSION(:):: zsumclasses ! 1D workspace 1574 !!------------------------------------------------------------- 1575 CALL wrk_alloc(nb_type_class , zsumclasses ) 1576 1577 zsumclasses(:)=0._wp 1578 zslope = sec%slopeSection 1579 1580 WRITE(numdct_NOOS,'(I4,a1,I2,a1,I2,a12,i3,a17,i3,a10,a25)') nyear,'.',nmonth,'.',nday,' Transect:',ksec-1,' No. of layers:',sec%nb_class-1,' Name: ',sec%name 1581 1582 DO jclass=1,MAX(1,sec%nb_class-1) 1583 zsumclasses(1:nb_type_class)=zsumclasses(1:nb_type_class)+sec%transport(1:nb_type_class,jclass) 1584 ENDDO 1585 1586 classe = 'total ' 1587 zbnd1 = 0._wp 1588 zbnd2 = 0._wp 1589 1590 IF ( zslope .gt. 0._wp .and. zslope .ne. 10000._wp ) THEN 1591 WRITE(numdct_NOOS,'(9e12.4E2)') -(zsumclasses( 1)+zsumclasses( 2)), -zsumclasses( 2),-zsumclasses( 1), & 1592 -(zsumclasses( 7)+zsumclasses( 8)), -zsumclasses( 8),-zsumclasses( 7), & 1593 -(zsumclasses( 9)+zsumclasses(10)), -zsumclasses(10),-zsumclasses( 9) 1594 ELSE 1595 WRITE(numdct_NOOS,'(9e12.4E2)') zsumclasses( 1)+zsumclasses( 2) , zsumclasses( 1), zsumclasses( 2), & 1596 zsumclasses( 7)+zsumclasses( 8) , zsumclasses( 7), zsumclasses( 8), & 1597 zsumclasses( 9)+zsumclasses(10) , zsumclasses( 9), zsumclasses(10) 1598 ENDIF 1599 1600 DO jclass=1,MAX(1,sec%nb_class-1) 1601 1602 classe = 'N ' 1603 zbnd1 = 0._wp 1604 zbnd2 = 0._wp 1605 1606 !insitu density classes transports 1607 IF( ( sec%zsigi(jclass) .NE. 99._wp ) .AND. & 1608 ( sec%zsigi(jclass+1) .NE. 99._wp ) )THEN 1609 classe = 'DI ' 1610 zbnd1 = sec%zsigi(jclass) 1611 zbnd2 = sec%zsigi(jclass+1) 1612 ENDIF 1613 !potential density classes transports 1614 IF( ( sec%zsigp(jclass) .NE. 99._wp ) .AND. & 1615 ( sec%zsigp(jclass+1) .NE. 99._wp ) )THEN 1616 classe = 'DP ' 1617 zbnd1 = sec%zsigp(jclass) 1618 zbnd2 = sec%zsigp(jclass+1) 1619 ENDIF 1620 !depth classes transports 1621 IF( ( sec%zlay(jclass) .NE. 99._wp ) .AND. & 1622 ( sec%zlay(jclass+1) .NE. 99._wp ) )THEN 1623 classe = 'Z ' 1624 zbnd1 = sec%zlay(jclass) 1625 zbnd2 = sec%zlay(jclass+1) 1626 ENDIF 1627 !salinity classes transports 1628 IF( ( sec%zsal(jclass) .NE. 99._wp ) .AND. & 1629 ( sec%zsal(jclass+1) .NE. 99._wp ) )THEN 1630 classe = 'S ' 1631 zbnd1 = sec%zsal(jclass) 1632 zbnd2 = sec%zsal(jclass+1) 1633 ENDIF 1634 !temperature classes transports 1635 IF( ( sec%ztem(jclass) .NE. 99._wp ) .AND. & 1636 ( sec%ztem(jclass+1) .NE. 99._wp ) ) THEN 1637 classe = 'T ' 1638 zbnd1 = sec%ztem(jclass) 1639 zbnd2 = sec%ztem(jclass+1) 1640 ENDIF 1641 1642 !write volume transport per class 1643 IF ( zslope .gt. 0._wp .and. zslope .ne. 10000._wp ) THEN 1644 WRITE(numdct_NOOS,'(9e12.4E2)') -(sec%transport( 1,jclass)+sec%transport( 2,jclass)),-sec%transport( 2,jclass),-sec%transport( 1,jclass), & 1645 -(sec%transport( 7,jclass)+sec%transport( 8,jclass)),-sec%transport( 8,jclass),-sec%transport( 7,jclass), & 1646 -(sec%transport( 9,jclass)+sec%transport(10,jclass)),-sec%transport(10,jclass),-sec%transport( 9,jclass) 1647 ELSE 1648 WRITE(numdct_NOOS,'(9e12.4E2)') sec%transport( 1,jclass)+sec%transport( 2,jclass) , sec%transport( 1,jclass), sec%transport( 2,jclass), & 1649 sec%transport( 7,jclass)+sec%transport( 8,jclass) , sec%transport( 7,jclass), sec%transport( 8,jclass), & 1650 sec%transport( 9,jclass)+sec%transport(10,jclass) , sec%transport( 9,jclass), sec%transport(10,jclass) 1651 ENDIF 1652 1653 ENDDO 1654 1655 CALL wrk_dealloc(nb_type_class , zsumclasses ) 1656 1657 END SUBROUTINE dia_dct_wri_NOOS 1658 1659 SUBROUTINE dia_dct_wri_NOOS_h(hr,ksec,sec) 1660 !!------------------------------------------------------------- 1661 !! As routine dia_dct_wri_NOOS but for hourly output files 1662 !! 1663 !! Write transport output in numdct using NOOS formatting 1664 !! 1665 !! Purpose: Write transports in ascii files 1666 !! 1667 !! Method: 1668 !! 1. Write volume transports in "volume_transport" 1669 !! Unit: Sv : area * Velocity / 1.e6 1670 !! 1671 !!------------------------------------------------------------- 1672 !!arguments 1673 INTEGER, INTENT(IN) :: hr ! hour 1674 TYPE(SECTION), INTENT(INOUT) :: sec ! section to write 1675 INTEGER ,INTENT(IN) :: ksec ! section number 1676 1677 !!local declarations 1678 INTEGER :: jclass,jhr ! Dummy loop 1679 CHARACTER(len=2) :: classe ! Classname 1680 REAL(wp) :: zbnd1,zbnd2 ! Class bounds 1681 REAL(wp) :: zslope ! section's slope coeff 1682 ! 1683 REAL(wp), POINTER, DIMENSION(:):: zsumclasses ! 1D workspace 1684 !!------------------------------------------------------------- 1685 1686 CALL wrk_alloc(nb_type_class , zsumclasses ) 1687 1688 zsumclasses(:)=0._wp 1689 zslope = sec%slopeSection 1690 1691 DO jclass=1,MAX(1,sec%nb_class-1) 1692 zsumclasses(1:nb_type_class)=zsumclasses(1:nb_type_class)+sec%transport_h(1:nb_type_class,jclass) 1693 ENDDO 1694 1695 !write volume transport per class 1696 IF ( zslope .gt. 0._wp .and. zslope .ne. 10000._wp ) THEN 1697 z_hr_output(ksec,hr,1)=-(zsumclasses(1)+zsumclasses(2)) 1698 ELSE 1699 z_hr_output(ksec,hr,1)= (zsumclasses(1)+zsumclasses(2)) 1700 ENDIF 1701 1702 DO jclass=1,MAX(1,sec%nb_class-1) 1703 IF ( zslope .gt. 0._wp .and. zslope .ne. 10000._wp ) THEN 1704 z_hr_output(ksec,hr,jclass+1)=-(sec%transport_h(1,jclass)+sec%transport_h(2,jclass)) 1705 ELSE 1706 z_hr_output(ksec,hr,jclass+1)= (sec%transport_h(1,jclass)+sec%transport_h(2,jclass)) 1707 ENDIF 1708 ENDDO 1709 1710 IF ( hr .eq. 48._wp ) THEN 1711 WRITE(numdct_NOOS_h,'(I4,a1,I2,a1,I2,a12,i3,a17,i3)') nyear,'.',nmonth,'.',nday,' Transect:',ksec-1,' No. of layers:',sec%nb_class-1 1712 DO jhr=25,48 1713 WRITE(numdct_NOOS_h,'(11F12.1)') z_hr_output(ksec,jhr,1), (z_hr_output(ksec,jhr,jclass+1), jclass=1,MAX(1,10) ) 1714 ENDDO 1715 ENDIF 1716 1717 CALL wrk_dealloc(nb_type_class , zsumclasses ) 1718 1719 END SUBROUTINE dia_dct_wri_NOOS_h 1720 1003 1721 SUBROUTINE dia_dct_wri(kt,ksec,sec) 1004 1722 !!------------------------------------------------------------- … … 1092 1810 WRITE(numdct_heat,119) ndastp,kt,ksec,sec%name,zslope, & 1093 1811 jclass,classe,zbnd1,zbnd2,& 1094 sec%transport( 3,jclass)*1.e-15,sec%transport(4,jclass)*1.e-15, &1095 ( sec%transport( 3,jclass)+sec%transport(4,jclass) )*1.e-151812 sec%transport(7,jclass)*1.e-15,sec%transport(8,jclass)*1.e-15, & 1813 ( sec%transport(7,jclass)+sec%transport(8,jclass) )*1.e-15 1096 1814 !write salt transport per class 1097 1815 WRITE(numdct_salt,119) ndastp,kt,ksec,sec%name,zslope, & 1098 1816 jclass,classe,zbnd1,zbnd2,& 1099 sec%transport( 5,jclass)*1.e-9,sec%transport(6,jclass)*1.e-9,&1100 (sec%transport( 5,jclass)+sec%transport(6,jclass))*1.e-91817 sec%transport(9,jclass)*1.e-9,sec%transport(10,jclass)*1.e-9,& 1818 (sec%transport(9,jclass)+sec%transport(10,jclass))*1.e-9 1101 1819 ENDIF 1102 1820 … … 1117 1835 WRITE(numdct_heat,119) ndastp,kt,ksec,sec%name,zslope, & 1118 1836 jclass,"total",zbnd1,zbnd2,& 1119 zsumclasses( 3)*1.e-15,zsumclasses(4)*1.e-15,&1120 (zsumclasses( 3)+zsumclasses(4) )*1.e-151837 zsumclasses(7)*1.e-15,zsumclasses(8)*1.e-15,& 1838 (zsumclasses(7)+zsumclasses(8) )*1.e-15 1121 1839 !write total salt transport 1122 1840 WRITE(numdct_salt,119) ndastp,kt,ksec,sec%name,zslope, & 1123 1841 jclass,"total",zbnd1,zbnd2,& 1124 zsumclasses( 5)*1.e-9,zsumclasses(6)*1.e-9,&1125 (zsumclasses( 5)+zsumclasses(6))*1.e-91842 zsumclasses(9)*1.e-9,zsumclasses(10)*1.e-9,& 1843 (zsumclasses(9)+zsumclasses(10))*1.e-9 1126 1844 ENDIF 1127 1845 … … 1131 1849 WRITE(numdct_vol,118) ndastp,kt,ksec,sec%name,zslope,& 1132 1850 jclass,"ice_vol",zbnd1,zbnd2,& 1133 sec%transport( 7,1),sec%transport(8,1),&1134 sec%transport( 7,1)+sec%transport(8,1)1851 sec%transport(11,1),sec%transport(12,1),& 1852 sec%transport(11,1)+sec%transport(12,1) 1135 1853 !write total ice surface transport 1136 1854 WRITE(numdct_vol,118) ndastp,kt,ksec,sec%name,zslope,& 1137 1855 jclass,"ice_surf",zbnd1,zbnd2,& 1138 sec%transport( 9,1),sec%transport(10,1), &1139 sec%transport( 9,1)+sec%transport(10,1)1856 sec%transport(13,1),sec%transport(14,1), & 1857 sec%transport(13,1)+sec%transport(14,1) 1140 1858 ENDIF 1141 1859 … … 1146 1864 END SUBROUTINE dia_dct_wri 1147 1865 1148 FUNCTION interp(ki, kj, kk, cd_point, ptab)1866 PURE FUNCTION interp(ki, kj, kk, cd_point,var) 1149 1867 !!---------------------------------------------------------------------- 1150 1868 !! … … 1208 1926 !*arguments 1209 1927 INTEGER, INTENT(IN) :: ki, kj, kk ! coordinate of point 1928 INTEGER, INTENT(IN) :: var ! which variable 1210 1929 CHARACTER(len=1), INTENT(IN) :: cd_point ! type of point (U, V) 1211 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(IN) :: ptab ! variable to compute at (ki, kj, kk )1212 1930 REAL(wp) :: interp ! interpolated variable 1213 1931 … … 1220 1938 !!---------------------------------------------------------------------- 1221 1939 1940 1941 1222 1942 IF( cd_point=='U' )THEN 1223 1943 ii1 = ki ; ij1 = kj … … 1250 1970 1251 1971 ! result 1252 interp = zmsk * ( zwgt2 * ptab(ii1,ij1,kk) + zwgt1 * ptab(ii1,ij1,kk) ) / ( zwgt2 + zwgt1 ) 1253 1972 SELECT CASE( var ) 1973 CASE(0) ; interp = zmsk * ( zwgt2 * tsn(ii1,ij1,kk,jp_tem) + zwgt1 * tsn(ii1,ij1,kk,jp_tem) ) / ( zwgt2 + zwgt1 ) 1974 CASE(1) ; interp = zmsk * ( zwgt2 * tsn(ii1,ij1,kk,jp_sal) + zwgt1 * tsn(ii1,ij1,kk,jp_sal) ) / ( zwgt2 + zwgt1 ) 1975 CASE(2) ; interp = zmsk * ( zwgt2 * rhop(ii1,ij1,kk) + zwgt1 * rhop(ii1,ij1,kk) ) / ( zwgt2 + zwgt1 ) 1976 CASE(3) ; interp = zmsk * ( zwgt2 * (rhd(ii1,ij1,kk)*rau0+rau0) + zwgt1 * (rhd(ii1,ij1,kk)*rau0+rau0) ) / ( zwgt2 + zwgt1 ) 1977 END SELECT 1254 1978 1255 1979 ELSE ! full step or partial step case … … 1273 1997 IF( ze3t >= 0. )THEN 1274 1998 ! zbis 1275 zbis = ptab(ii2,ij2,kk) + zwgt1 * ( ptab(ii2,ij2,kk-1) - ptab(ii2,ij2,kk) ) 1999 SELECT CASE( var ) 2000 CASE(0) 2001 zbis = tsn(ii2,ij2,kk,jp_tem) + zwgt1 * (tsn(ii2,ij2,kk-1,jp_tem)-tsn(ii2,ij2,kk,jp_tem) ) 2002 interp = zmsk * ( zet2 * tsn(ii1,ij1,kk,jp_tem) + zet1 * zbis )/( zet1 + zet2 ) 2003 CASE(1) 2004 zbis = tsn(ii2,ij2,kk,jp_sal) + zwgt1 * (tsn(ii2,ij2,kk-1,jp_sal)-tsn(ii2,ij2,kk,jp_sal) ) 2005 interp = zmsk * ( zet2 * tsn(ii1,ij1,kk,jp_sal) + zet1 * zbis )/( zet1 + zet2 ) 2006 CASE(2) 2007 zbis = rhop(ii2,ij2,kk) + zwgt1 * (rhop(ii2,ij2,kk-1)-rhop(ii2,ij2,kk) ) 2008 interp = zmsk * ( zet2 * rhop(ii1,ij1,kk) + zet1 * zbis )/( zet1 + zet2 ) 2009 CASE(3) 2010 zbis = (rhd(ii2,ij2,kk)*rau0+rau0) + zwgt1 * ( (rhd(ii2,ij2,kk-1)*rau0+rau0)-(rhd(ii2,ij2,kk)*rau0+rau0) ) 2011 interp = zmsk * ( zet2 * (rhd(ii1,ij1,kk)*rau0+rau0) + zet1 * zbis )/( zet1 + zet2 ) 2012 END SELECT 1276 2013 ! result 1277 interp = zmsk * ( zet2 * ptab(ii1,ij1,kk) + zet1 * zbis )/( zet1 + zet2 )1278 2014 ELSE 1279 2015 ! zbis 1280 zbis = ptab(ii1,ij1,kk) + zwgt2 * ( ptab(ii1,ij1,kk-1) - ptab(ii1,ij2,kk) ) 1281 ! result 1282 interp = zmsk * ( zet2 * zbis + zet1 * ptab(ii2,ij2,kk) )/( zet1 + zet2 ) 2016 SELECT CASE( var ) 2017 CASE(0) 2018 zbis = tsn(ii1,ij1,kk,jp_tem) + zwgt2 * ( tsn(ii1,ij1,kk-1,jp_tem) - tsn(ii1,ij2,kk,jp_tem) ) 2019 interp = zmsk * ( zet2 * zbis + zet1 * tsn(ii2,ij2,kk,jp_tem) )/( zet1 + zet2 ) 2020 CASE(1) 2021 zbis = tsn(ii1,ij1,kk,jp_sal) + zwgt2 * ( tsn(ii1,ij1,kk-1,jp_sal) - tsn(ii1,ij2,kk,jp_sal) ) 2022 interp = zmsk * ( zet2 * zbis + zet1 * tsn(ii2,ij2,kk,jp_sal) )/( zet1 + zet2 ) 2023 CASE(2) 2024 zbis = rhop(ii1,ij1,kk) + zwgt2 * ( rhop(ii1,ij1,kk-1) - rhop(ii1,ij2,kk) ) 2025 interp = zmsk * ( zet2 * zbis + zet1 * rhop(ii2,ij2,kk) )/( zet1 + zet2 ) 2026 CASE(3) 2027 zbis = (rhd(ii1,ij1,kk)*rau0+rau0) + zwgt2 * ( (rhd(ii1,ij1,kk-1)*rau0+rau0) - (rhd(ii1,ij2,kk)*rau0+rau0) ) 2028 interp = zmsk * ( zet2 * zbis + zet1 * (rhd(ii2,ij2,kk)*rau0+rau0) )/( zet1 + zet2 ) 2029 END SELECT 1283 2030 ENDIF 1284 2031 1285 2032 ELSE 1286 interp = zmsk * ( zet2 * ptab(ii1,ij1,kk) + zet1 * ptab(ii2,ij2,kk) )/( zet1 + zet2 ) 2033 SELECT CASE( var ) 2034 CASE(0) 2035 interp = zmsk * ( zet2 * tsn(ii1,ij1,kk,jp_tem) + zet1 * tsn(ii2,ij2,kk,jp_tem) )/( zet1 + zet2 ) 2036 CASE(1) 2037 interp = zmsk * ( zet2 * tsn(ii1,ij1,kk,jp_sal) + zet1 * tsn(ii2,ij2,kk,jp_sal) )/( zet1 + zet2 ) 2038 CASE(2) 2039 interp = zmsk * ( zet2 * rhop(ii1,ij1,kk) + zet1 * rhop(ii2,ij2,kk) )/( zet1 + zet2 ) 2040 CASE(3) 2041 interp = zmsk * ( zet2 * (rhd(ii1,ij1,kk)*rau0+rau0) + zet1 * (rhd(ii2,ij2,kk)*rau0+rau0) )/( zet1 + zet2 ) 2042 END SELECT 1287 2043 ENDIF 1288 2044 1289 2045 ENDIF 1290 1291 2046 1292 2047 END FUNCTION interp -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
r8059 r10237 135 135 !!---------------------------------------------------------------------- 136 136 USE ioipsl 137 NAMELIST/namrun/ ln_NOOS 137 138 NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list, & 138 139 & nn_no , cn_exp , cn_ocerst_in, cn_ocerst_out, ln_rstdate, ln_rstart , nn_rstctl, & … … 193 194 WRITE(numout,*) ' overwrite an existing file ln_clobber = ', ln_clobber 194 195 WRITE(numout,*) ' NetCDF chunksize (bytes) nn_chunksz = ', nn_chunksz 196 WRITE(numout,*) ' NOOS transect diagnostics ln_NOOS = ', ln_NOOS 195 197 ENDIF 196 198 -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/IOM/in_out_manager.F90
r8059 r10237 49 49 LOGICAL :: ln_clobber !: clobber (overwrite) an existing file 50 50 INTEGER :: nn_chunksz !: chunksize (bytes) for NetCDF file (works only with iom_nf90 routines) 51 LOGICAL :: ln_NOOS = .FALSE. !: NOOS transects diagnostics 51 52 #if defined key_netcdf4 52 53 !!---------------------------------------------------------------------- … … 134 135 INTEGER :: numdct_heat = -1 !: logical unit for heat transports output 135 136 INTEGER :: numdct_salt = -1 !: logical unit for salt transports output 137 INTEGER :: numdct_NOOS = -1 !: logical unit for NOOS transports output 138 INTEGER :: numdct_NOOS_h = -1 !: logical unit for NOOS hourly transports output 136 139 INTEGER :: numfl = -1 !: logical unit for floats ascii output 137 140 INTEGER :: numflo = -1 !: logical unit for floats ascii output 141 142 INTEGER :: numdct_reg_bin = -1 !: logical unit for NOOS transports output 143 INTEGER :: numdct_reg_txt = -1 !: logical unit for NOOS hourly transports output 138 144 139 145 !!----------------------------------------------------------------------
Note: See TracChangeset
for help on using the changeset viewer.