Changeset 6140 for trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90
- Timestamp:
- 2015-12-21T12:35:23+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90
r5505 r6140 1 1 MODULE diadct 2 !!===================================================================== 3 !! *** MODULE diadct *** 4 !! Ocean diagnostics: Compute the transport trough a sec. 5 !!=============================================================== 6 !! History : 7 !! 8 !! original : 02/99 (Y Drillet) 9 !! addition : 10/01 (Y Drillet, R Bourdalle Badie) 10 !! : 10/05 (M Laborie) F90 11 !! addition : 04/07 (G Garric) Ice sections 12 !! bugfix : 04/07 (C Bricaud) test on sec%nb_point 13 !! initialisation of ztransp1,ztransp2,... 14 !! nemo_v_3_4: 09/2011 (C Bricaud) 15 !! 16 !! 17 !!---------------------------------------------------------------------- 2 !!====================================================================== 3 !! *** MODULE diadct *** 4 !! Ocean diagnostics: Compute the transport trough a sec. 5 !!====================================================================== 6 !! History : OPA ! 02/1999 (Y Drillet) original code 7 !! ! 10/2001 (Y Drillet, R Bourdalle Badie) 8 !! NEMO 1.0 ! 10/2005 (M Laborie) F90 9 !! 3.0 ! 04/2007 (G Garric) Ice sections 10 !! - ! 04/2007 (C Bricaud) test on sec%nb_point, initialisation of ztransp1,ztransp2,... 11 !! 3.4 ! 09/2011 (C Bricaud) 12 !!---------------------------------------------------------------------- 18 13 #if defined key_diadct 19 !!---------------------------------------------------------------------- 20 !! 'key_diadct' : 21 !!---------------------------------------------------------------------- 22 !!---------------------------------------------------------------------- 23 !! dia_dct : Compute the transport through a sec. 24 !! dia_dct_init : Read namelist. 25 !! readsec : Read sections description and pathway 26 !! removepoints : Remove points which are common to 2 procs 27 !! transport : Compute transport for each sections 28 !! dia_dct_wri : Write tranports results in ascii files 29 !! interp : Compute temperature/salinity/density at U-point or V-point 30 !! 31 !!---------------------------------------------------------------------- 32 !! * Modules used 33 USE oce ! ocean dynamics and tracers 34 USE dom_oce ! ocean space and time domain 35 USE phycst ! physical constants 36 USE in_out_manager ! I/O manager 37 USE daymod ! calendar 38 USE dianam ! build name of file 39 USE lib_mpp ! distributed memory computing library 14 !!---------------------------------------------------------------------- 15 !! 'key_diadct' : 16 !!---------------------------------------------------------------------- 17 !!---------------------------------------------------------------------- 18 !! dia_dct : Compute the transport through a sec. 19 !! dia_dct_init : Read namelist. 20 !! readsec : Read sections description and pathway 21 !! removepoints : Remove points which are common to 2 procs 22 !! transport : Compute transport for each sections 23 !! dia_dct_wri : Write tranports results in ascii files 24 !! interp : Compute temperature/salinity/density at U-point or V-point 25 !! 26 !!---------------------------------------------------------------------- 27 USE oce ! ocean dynamics and tracers 28 USE dom_oce ! ocean space and time domain 29 USE phycst ! physical constants 30 USE in_out_manager ! I/O manager 31 USE daymod ! calendar 32 USE dianam ! build name of file 33 USE lib_mpp ! distributed memory computing library 40 34 #if defined key_lim2 41 USE ice_235 USE ice_2 42 36 #endif 43 37 #if defined key_lim3 44 USE ice38 USE ice 45 39 #endif 46 USE domvvl 47 USE timing ! preformance summary 48 USE wrk_nemo ! working arrays 49 50 IMPLICIT NONE 51 PRIVATE 52 53 !! * Routine accessibility 54 PUBLIC dia_dct ! routine called by step.F90 55 PUBLIC dia_dct_init ! routine called by opa.F90 56 PUBLIC diadct_alloc ! routine called by nemo_init in nemogcm.F90 57 PRIVATE readsec 58 PRIVATE removepoints 59 PRIVATE transport 60 PRIVATE dia_dct_wri 61 62 #include "domzgr_substitute.h90" 63 64 !! * Shared module variables 65 LOGICAL, PUBLIC, PARAMETER :: lk_diadct = .TRUE. !: model-data diagnostics flag 66 67 !! * 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 40 USE domvvl 41 USE timing ! preformance summary 42 USE wrk_nemo ! working arrays 43 44 IMPLICIT NONE 45 PRIVATE 46 47 PUBLIC dia_dct ! routine called by step.F90 48 PUBLIC dia_dct_init ! routine called by opa.F90 49 PUBLIC diadct_alloc ! routine called by nemo_init in nemogcm.F90 50 PRIVATE readsec 51 PRIVATE removepoints 52 PRIVATE transport 53 PRIVATE dia_dct_wri 54 55 LOGICAL, PUBLIC, PARAMETER :: lk_diadct = .TRUE. !: model-data diagnostics flag 56 57 INTEGER :: nn_dct ! Frequency of computation 58 INTEGER :: nn_dctwri ! Frequency of output 59 INTEGER :: nn_secdebug ! Number of the section to debug 71 60 72 INTEGER, PARAMETER :: nb_class_max = 1073 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 = 278 INTEGER :: nb_sec79 80 TYPE POINT_SECTION81 INTEGER :: I,J82 END TYPE POINT_SECTION83 84 TYPE COORD_SECTION85 REAL(wp) :: lon,lat86 END TYPE COORD_SECTION87 88 TYPE SECTION89 CHARACTER(len=60) :: name ! name of the sec90 LOGICAL :: llstrpond ! true if you want the computation of salt and61 INTEGER, PARAMETER :: nb_class_max = 10 62 INTEGER, PARAMETER :: nb_sec_max = 150 63 INTEGER, PARAMETER :: nb_point_max = 2000 64 INTEGER, PARAMETER :: nb_type_class = 10 65 INTEGER, PARAMETER :: nb_3d_vars = 3 66 INTEGER, PARAMETER :: nb_2d_vars = 2 67 INTEGER :: nb_sec 68 69 TYPE POINT_SECTION 70 INTEGER :: I,J 71 END TYPE POINT_SECTION 72 73 TYPE COORD_SECTION 74 REAL(wp) :: lon,lat 75 END TYPE COORD_SECTION 76 77 TYPE SECTION 78 CHARACTER(len=60) :: name ! name of the sec 79 LOGICAL :: llstrpond ! true if you want the computation of salt and 91 80 ! heat transports 92 LOGICAL :: ll_ice_section ! ice surface and ice volume computation 93 LOGICAL :: ll_date_line ! = T if the section crosses the date-line 94 TYPE(COORD_SECTION), DIMENSION(2) :: coordSec ! longitude and latitude of the extremities of the sec 95 INTEGER :: nb_class ! number of boundaries for density classes 96 INTEGER, DIMENSION(nb_point_max) :: direction ! vector direction of the point in the section 97 CHARACTER(len=40),DIMENSION(nb_class_max) :: classname ! characteristics of the class 98 REAL(wp), DIMENSION(nb_class_max) :: zsigi ,&! in-situ density classes (99 if you don't want) 99 zsigp ,&! potential density classes (99 if you don't want) 100 zsal ,&! salinity classes (99 if you don't want) 101 ztem ,&! temperature classes(99 if you don't want) 102 zlay ! level classes (99 if you don't want) 103 REAL(wp), DIMENSION(nb_type_class,nb_class_max) :: transport ! transport output 104 REAL(wp) :: slopeSection ! slope of the section 105 INTEGER :: nb_point ! number of points in the section 106 TYPE(POINT_SECTION),DIMENSION(nb_point_max) :: listPoint ! list of points in the sections 107 END TYPE SECTION 108 109 TYPE(SECTION),DIMENSION(nb_sec_max) :: secs ! Array of sections 110 111 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: transports_3d 112 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: transports_2d 113 81 LOGICAL :: ll_ice_section ! ice surface and ice volume computation 82 LOGICAL :: ll_date_line ! = T if the section crosses the date-line 83 TYPE(COORD_SECTION), DIMENSION(2) :: coordSec ! longitude and latitude of the extremities of the sec 84 INTEGER :: nb_class ! number of boundaries for density classes 85 INTEGER, DIMENSION(nb_point_max) :: direction ! vector direction of the point in the section 86 CHARACTER(len=40),DIMENSION(nb_class_max) :: classname ! characteristics of the class 87 REAL(wp), DIMENSION(nb_class_max) :: zsigi ,&! in-situ density classes (99 if you don't want) 88 zsigp ,&! potential density classes (99 if you don't want) 89 zsal ,&! salinity classes (99 if you don't want) 90 ztem ,&! temperature classes(99 if you don't want) 91 zlay ! level classes (99 if you don't want) 92 REAL(wp), DIMENSION(nb_type_class,nb_class_max) :: transport ! transport output 93 REAL(wp) :: slopeSection ! slope of the section 94 INTEGER :: nb_point ! number of points in the section 95 TYPE(POINT_SECTION),DIMENSION(nb_point_max) :: listPoint ! list of points in the sections 96 END TYPE SECTION 97 98 TYPE(SECTION),DIMENSION(nb_sec_max) :: secs ! Array of sections 99 100 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: transports_3d 101 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: transports_2d 102 103 !!---------------------------------------------------------------------- 104 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 114 105 !! $Id$ 106 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 107 !!---------------------------------------------------------------------- 115 108 CONTAINS 116 117 109 118 110 INTEGER FUNCTION diadct_alloc() … … 130 122 131 123 END FUNCTION diadct_alloc 124 132 125 133 126 SUBROUTINE dia_dct_init … … 208 201 !! Reinitialise all relevant arrays to zero 209 202 !!--------------------------------------------------------------------- 210 !! * Arguments 211 INTEGER,INTENT(IN) ::kt 212 213 !! * Local variables 203 INTEGER,INTENT(in) ::kt 204 ! 214 205 INTEGER :: jsec, &! loop on sections 215 206 itotal ! nb_sec_max*nb_type_class*nb_class_max … … 220 211 REAL(wp), POINTER, DIMENSION(:) :: zwork ! " 221 212 REAL(wp), POINTER, DIMENSION(:,:,:):: zsum ! " 222 223 213 !!--------------------------------------------------------------------- 214 ! 224 215 IF( nn_timing == 1 ) CALL timing_start('dia_dct') 225 216 … … 619 610 zumid_ice, zvmid_ice, &!U/V ice velocity 620 611 zTnorm !transport of velocity through one cell's sides 621 REAL(wp) :: ztn, zsn, zrhoi, zrhop, zsshn, z fsdep !temperature/salinity/potential density/ssh/depth at u/v point612 REAL(wp) :: ztn, zsn, zrhoi, zrhop, zsshn, zdep !temperature/salinity/potential density/ssh/depth at u/v point 622 613 623 614 TYPE(POINT_SECTION) :: k 624 !!--------------------------------------------------------625 626 IF( ld_debug )WRITE(numout,*)' Compute transport'627 628 !---------------------------!629 ! COMPUTE TRANSPORT !630 !---------------------------!631 IF(sec%nb_point .NE. 0)THEN632 633 !----------------------------------------------------------------------------------------------------634 !Compute sign for velocities:635 !636 !convention:637 ! non horizontal section: direction + is toward left hand of section638 ! horizontal section: direction + is toward north of section639 !640 !641 ! slopeSection < 0 slopeSection > 0 slopeSection=inf slopeSection=0642 ! ---------------- ----------------- --------------- --------------643 !644 ! isgnv=1 direction +645 ! ______ _____ ______646 ! | //| | | direction +647 ! | isgnu=1 // | |isgnu=1 |isgnu=1 /|\648 ! |_______ // ______| \\ | ---\ |649 ! | | isgnv=-1 \\ | | ---/ direction + ____________650 ! | | __\\| |651 ! | | direction + | isgnv=1652 !653 !----------------------------------------------------------------------------------------------------654 isgnu = 1655 IF( sec%slopeSection .GT. 0 ) THEN ; isgnv = -1656 ELSE ; isgnv = 1657 ENDIF658 IF( sec%slopeSection .GE. 9999. ) isgnv = 1659 660 IF( ld_debug )write(numout,*)"sec%slopeSection isgnu isgnv ",sec%slopeSection,isgnu,isgnv661 662 !--------------------------------------!663 ! LOOP ON THE SEGMENT BETWEEN 2 NODES !664 !--------------------------------------!665 DO jseg=1,MAX(sec%nb_point-1,0)615 !!-------------------------------------------------------- 616 ! 617 IF( ld_debug )WRITE(numout,*)' Compute transport' 618 619 !---------------------------! 620 ! COMPUTE TRANSPORT ! 621 !---------------------------! 622 IF(sec%nb_point .NE. 0)THEN 623 624 !---------------------------------------------------------------------------------------------------- 625 !Compute sign for velocities: 626 ! 627 !convention: 628 ! non horizontal section: direction + is toward left hand of section 629 ! horizontal section: direction + is toward north of section 630 ! 631 ! 632 ! slopeSection < 0 slopeSection > 0 slopeSection=inf slopeSection=0 633 ! ---------------- ----------------- --------------- -------------- 634 ! 635 ! isgnv=1 direction + 636 ! ______ _____ ______ 637 ! | //| | | direction + 638 ! | isgnu=1 // | |isgnu=1 |isgnu=1 /|\ 639 ! |_______ // ______| \\ | ---\ | 640 ! | | isgnv=-1 \\ | | ---/ direction + ____________ 641 ! | | __\\| | 642 ! | | direction + | isgnv=1 643 ! 644 !---------------------------------------------------------------------------------------------------- 645 isgnu = 1 646 IF( sec%slopeSection .GT. 0 ) THEN ; isgnv = -1 647 ELSE ; isgnv = 1 648 ENDIF 649 IF( sec%slopeSection .GE. 9999. ) isgnv = 1 650 651 IF( ld_debug )write(numout,*)"sec%slopeSection isgnu isgnv ",sec%slopeSection,isgnu,isgnv 652 653 !--------------------------------------! 654 ! LOOP ON THE SEGMENT BETWEEN 2 NODES ! 655 !--------------------------------------! 656 DO jseg=1,MAX(sec%nb_point-1,0) 666 657 667 !------------------------------------------------------------------------------------------- 668 ! Select the appropriate coordinate for computing the velocity of the segment 669 ! 670 ! CASE(0) Case (2) 671 ! ------- -------- 672 ! listPoint(jseg) listPoint(jseg+1) listPoint(jseg) F(i,j) 673 ! F(i,j)----------V(i+1,j)-------F(i+1,j) | 674 ! | 675 ! | 676 ! | 677 ! Case (3) U(i,j) 678 ! -------- | 679 ! | 680 ! listPoint(jseg+1) F(i,j+1) | 681 ! | | 682 ! | | 683 ! | listPoint(jseg+1) F(i,j-1) 684 ! | 685 ! | 686 ! U(i,j+1) 687 ! | Case(1) 688 ! | ------ 689 ! | 690 ! | listPoint(jseg+1) listPoint(jseg) 691 ! | F(i-1,j)-----------V(i,j) -------f(jseg) 692 ! listPoint(jseg) F(i,j) 693 ! 694 !------------------------------------------------------------------------------------------- 695 696 SELECT CASE( sec%direction(jseg) ) 697 CASE(0) ; k = sec%listPoint(jseg) 698 CASE(1) ; k = POINT_SECTION(sec%listPoint(jseg)%I+1,sec%listPoint(jseg)%J) 699 CASE(2) ; k = sec%listPoint(jseg) 700 CASE(3) ; k = POINT_SECTION(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J+1) 701 END SELECT 702 703 !---------------------------| 704 ! LOOP ON THE LEVEL | 705 !---------------------------| 706 !Sum of the transport on the vertical 707 DO jk=1,mbathy(k%I,k%J) 708 709 ! compute temperature, salinity, insitu & potential density, ssh and depth at U/V point 710 SELECT CASE( sec%direction(jseg) ) 711 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) 716 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) * vmask(k%I,k%J,1) 717 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) 722 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) * umask(k%I,k%J,1) 723 END SELECT 724 725 zfsdep= fsdept(k%I,k%J,jk) 658 !------------------------------------------------------------------------------------------- 659 ! Select the appropriate coordinate for computing the velocity of the segment 660 ! 661 ! CASE(0) Case (2) 662 ! ------- -------- 663 ! listPoint(jseg) listPoint(jseg+1) listPoint(jseg) F(i,j) 664 ! F(i,j)----------V(i+1,j)-------F(i+1,j) | 665 ! | 666 ! | 667 ! | 668 ! Case (3) U(i,j) 669 ! -------- | 670 ! | 671 ! listPoint(jseg+1) F(i,j+1) | 672 ! | | 673 ! | | 674 ! | listPoint(jseg+1) F(i,j-1) 675 ! | 676 ! | 677 ! U(i,j+1) 678 ! | Case(1) 679 ! | ------ 680 ! | 681 ! | listPoint(jseg+1) listPoint(jseg) 682 ! | F(i-1,j)-----------V(i,j) -------f(jseg) 683 ! listPoint(jseg) F(i,j) 684 ! 685 !------------------------------------------------------------------------------------------- 686 687 SELECT CASE( sec%direction(jseg) ) 688 CASE(0) ; k = sec%listPoint(jseg) 689 CASE(1) ; k = POINT_SECTION(sec%listPoint(jseg)%I+1,sec%listPoint(jseg)%J) 690 CASE(2) ; k = sec%listPoint(jseg) 691 CASE(3) ; k = POINT_SECTION(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J+1) 692 END SELECT 693 694 !---------------------------| 695 ! LOOP ON THE LEVEL | 696 !---------------------------| 697 DO jk = 1, mbathy(k%I,k%J) !Sum of the transport on the vertical 698 ! ! compute temperature, salinity, insitu & potential density, ssh and depth at U/V point 699 SELECT CASE( sec%direction(jseg) ) 700 CASE(0,1) 701 ztn = interp(k%I,k%J,jk,'V',tsn(:,:,:,jp_tem) ) 702 zsn = interp(k%I,k%J,jk,'V',tsn(:,:,:,jp_sal) ) 703 zrhop = interp(k%I,k%J,jk,'V',rhop) 704 zrhoi = interp(k%I,k%J,jk,'V',rhd*rau0+rau0) 705 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I,k%J+1) ) * vmask(k%I,k%J,1) 706 CASE(2,3) 707 ztn = interp(k%I,k%J,jk,'U',tsn(:,:,:,jp_tem) ) 708 zsn = interp(k%I,k%J,jk,'U',tsn(:,:,:,jp_sal) ) 709 zrhop = interp(k%I,k%J,jk,'U',rhop) 710 zrhoi = interp(k%I,k%J,jk,'U',rhd*rau0+rau0) 711 zsshn = 0.5*( sshn(k%I,k%J) + sshn(k%I+1,k%J) ) * umask(k%I,k%J,1) 712 END SELECT 713 ! 714 zdep= gdept_n(k%I,k%J,jk) 726 715 727 !compute velocity with the correct direction728 SELECT CASE( sec%direction(jseg) )729 CASE(0,1)730 zumid=0.731 zvmid=isgnv*vn(k%I,k%J,jk)*vmask(k%I,k%J,jk)732 CASE(2,3)733 zumid=isgnu*un(k%I,k%J,jk)*umask(k%I,k%J,jk)734 zvmid=0.735 END SELECT736 737 !zTnorm=transport through one cell;738 !velocity* cell's length * cell's thickness739 zTnorm=zumid*e2u(k%I,k%J)* fse3u(k%I,k%J,jk)+ &740 zvmid*e1v(k%I,k%J)* fse3v(k%I,k%J,jk) 741 742 #if ! defined key_vvl 743 !add transport due to free surface744 IF( jk==1 )THEN745 zTnorm = zTnorm + zumid* e2u(k%I,k%J) * zsshn * umask(k%I,k%J,jk) + &746 zvmid* e1v(k%I,k%J) * zsshn * vmask(k%I,k%J,jk)747 ENDIF748 #endif 716 SELECT CASE( sec%direction(jseg) ) !compute velocity with the correct direction 717 CASE(0,1) 718 zumid=0._wp 719 zvmid=isgnv*vn(k%I,k%J,jk)*vmask(k%I,k%J,jk) 720 CASE(2,3) 721 zumid=isgnu*un(k%I,k%J,jk)*umask(k%I,k%J,jk) 722 zvmid=0._wp 723 END SELECT 724 725 !zTnorm=transport through one cell; 726 !velocity* cell's length * cell's thickness 727 zTnorm = zumid*e2u(k%I,k%J) * e3u_n(k%I,k%J,jk) & 728 & + zvmid*e1v(k%I,k%J) * e3v_n(k%I,k%J,jk) 729 730 !!gm THIS is WRONG no transport due to ssh in linear free surface case !!!!! 731 IF( ln_linssh ) THEN !add transport due to free surface 732 IF( jk==1 ) THEN 733 zTnorm = zTnorm + zumid* e2u(k%I,k%J) * zsshn * umask(k%I,k%J,jk) & 734 & + zvmid* e1v(k%I,k%J) * zsshn * vmask(k%I,k%J,jk) 735 ENDIF 736 ENDIF 737 !!gm end 749 738 !COMPUTE TRANSPORT 750 739 751 740 transports_3d(1,jsec,jseg,jk) = transports_3d(1,jsec,jseg,jk) + zTnorm 752 741 753 IF 742 IF( sec%llstrpond ) THEN 754 743 transports_3d(2,jsec,jseg,jk) = transports_3d(2,jsec,jseg,jk) + zTnorm * ztn * zrhop * rcp 755 744 transports_3d(3,jsec,jseg,jk) = transports_3d(3,jsec,jseg,jk) + zTnorm * zsn * zrhop * 0.001 756 745 ENDIF 757 746 758 END DO !end of loop on the level747 END DO !end of loop on the level 759 748 760 749 #if defined key_lim2 || defined key_lim3 … … 797 786 transports_2d(2,jsec,jseg) = transports_2d(2,jsec,jseg) + (zTnorm)* & 798 787 a_i(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl) 799 END DO788 END DO 800 789 #endif 801 790 … … 803 792 #endif 804 793 805 END DO !end of loop on the segment794 END DO !end of loop on the segment 806 795 807 796 ENDIF !end of sec%nb_point =0 case 808 797 ! 809 798 END SUBROUTINE transport 810 799 800 811 801 SUBROUTINE dia_dct_sum(sec,jsec) 812 802 !!------------------------------------------------------------- … … 828 818 !! 829 819 !!------------------------------------------------------------- 830 !! * arguments831 820 TYPE(SECTION),INTENT(INOUT) :: sec 832 821 INTEGER ,INTENT(IN) :: jsec ! numeric identifier of section … … 834 823 TYPE(POINT_SECTION) :: k 835 824 INTEGER :: jk,jseg,jclass ! dummy variables for looping on level/segment/classes 836 REAL(wp) :: ztn, zsn, zrhoi, zrhop, zsshn, z fsdep ! temperature/salinity/ssh/potential density /depth at u/v point825 REAL(wp) :: ztn, zsn, zrhoi, zrhop, zsshn, zdep ! temperature/salinity/ssh/potential density /depth at u/v point 837 826 !!------------------------------------------------------------- 838 827 … … 903 892 END SELECT 904 893 905 z fsdep= fsdept(k%I,k%J,jk)894 zdep= gdept_n(k%I,k%J,jk) 906 895 907 896 !------------------------------- … … 932 921 ( sec%ztem(jclass) .EQ.99.)) .AND. & 933 922 934 ((( z fsdep .GE. sec%zlay(jclass)) .AND. &935 ( z fsdep .LE. sec%zlay(jclass+1))) .OR. &923 ((( zdep .GE. sec%zlay(jclass)) .AND. & 924 ( zdep .LE. sec%zlay(jclass+1))) .OR. & 936 925 ( sec%zlay(jclass) .EQ. 99. )) & 937 926 )) THEN … … 991 980 #endif 992 981 993 END DO !end of loop on the segment982 END DO !end of loop on the segment 994 983 995 984 ELSE !if sec%nb_point =0 … … 1000 989 1001 990 END SUBROUTINE dia_dct_sum 1002 991 992 1003 993 SUBROUTINE dia_dct_wri(kt,ksec,sec) 1004 994 !!------------------------------------------------------------- … … 1138 1128 sec%transport(9,1),sec%transport(10,1), & 1139 1129 sec%transport(9,1)+sec%transport(10,1) 1140 ENDIF1130 ENDIF 1141 1131 1142 118 FORMAT(I8,1X,I8,1X,I4,1X,A30,1X,f9.2,1X,I4,3X,A8,1X,2F12.4,5X,3F12.4) 1143 119 FORMAT(I8,1X,I8,1X,I4,1X,A30,1X,f9.2,1X,I4,3X,A8,1X,2F12.4,5X,3E15.6) 1144 1145 CALL wrk_dealloc(nb_type_class , zsumclasses ) 1146 END SUBROUTINE dia_dct_wri 1147 1148 FUNCTION interp(ki, kj, kk, cd_point, ptab) 1132 118 FORMAT(I8,1X,I8,1X,I4,1X,A30,1X,f9.2,1X,I4,3X,A8,1X,2F12.4,5X,3F12.4) 1133 119 FORMAT(I8,1X,I8,1X,I4,1X,A30,1X,f9.2,1X,I4,3X,A8,1X,2F12.4,5X,3E15.6) 1134 1135 CALL wrk_dealloc(nb_type_class , zsumclasses ) 1136 ! 1137 END SUBROUTINE dia_dct_wri 1138 1139 1140 FUNCTION interp(ki, kj, kk, cd_point, ptab) 1149 1141 !!---------------------------------------------------------------------- 1150 1142 !! … … 1214 1206 !*local declations 1215 1207 INTEGER :: ii1, ij1, ii2, ij2 ! local integer 1216 REAL(wp):: ze3t, z fse3, zwgt1, zwgt2, zbis, zdepu ! local real1208 REAL(wp):: ze3t, ze3, zwgt1, zwgt2, zbis, zdepu ! local real 1217 1209 REAL(wp):: zet1, zet2 ! weight for interpolation 1218 1210 REAL(wp):: zdep1,zdep2 ! differences of depth … … 1241 1233 IF( ln_sco )THEN ! s-coordinate case 1242 1234 1243 zdepu = ( fsdept(ii1,ij1,kk) + fsdept(ii2,ij2,kk) ) /21244 zdep1 = fsdept(ii1,ij1,kk) - zdepu1245 zdep2 = fsdept(ii2,ij2,kk) - zdepu1235 zdepu = ( gdept_n(ii1,ij1,kk) + gdept_n(ii2,ij2,kk) ) * 0.5_wp 1236 zdep1 = gdept_n(ii1,ij1,kk) - zdepu 1237 zdep2 = gdept_n(ii2,ij2,kk) - zdepu 1246 1238 1247 1239 ! weights … … 1255 1247 ELSE ! full step or partial step case 1256 1248 1257 #if defined key_vvl 1258 1259 ze3t = fse3t_n(ii2,ij2,kk) - fse3t_n(ii1,ij1,kk) 1260 zwgt1 = ( fse3w_n(ii2,ij2,kk) - fse3w_n(ii1,ij1,kk) ) / fse3w_n(ii2,ij2,kk) 1261 zwgt2 = ( fse3w_n(ii1,ij1,kk) - fse3w_n(ii2,ij2,kk) ) / fse3w_n(ii1,ij1,kk) 1262 1263 #else 1264 1265 ze3t = fse3t(ii2,ij2,kk) - fse3t(ii1,ij1,kk) 1266 zwgt1 = ( fse3w(ii2,ij2,kk) - fse3w(ii1,ij1,kk) ) / fse3w(ii2,ij2,kk) 1267 zwgt2 = ( fse3w(ii1,ij1,kk) - fse3w(ii2,ij2,kk) ) / fse3w(ii1,ij1,kk) 1268 1269 #endif 1249 ze3t = e3t_n(ii2,ij2,kk) - e3t_n(ii1,ij1,kk) 1250 zwgt1 = ( e3w_n(ii2,ij2,kk) - e3w_n(ii1,ij1,kk) ) / e3w_n(ii2,ij2,kk) 1251 zwgt2 = ( e3w_n(ii1,ij1,kk) - e3w_n(ii2,ij2,kk) ) / e3w_n(ii1,ij1,kk) 1270 1252 1271 1253 IF(kk .NE. 1)THEN … … 1288 1270 1289 1271 ENDIF 1290 1291 1292 END FUNCTION interp 1272 ! 1273 END FUNCTION interp 1293 1274 1294 1275 #else … … 1311 1292 #endif 1312 1293 1294 !!====================================================================== 1313 1295 END MODULE diadct
Note: See TracChangeset
for help on using the changeset viewer.