New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 13463 for NEMO/branches/2019/dev_r11351_fldread_with_XIOS/src/OCE/SBC/sbcdcy.F90 – NEMO

Ignore:
Timestamp:
2020-09-14T17:40:34+02:00 (4 years ago)
Author:
andmirek
Message:

Ticket #2195:update to trunk 13461

Location:
NEMO/branches/2019/dev_r11351_fldread_with_XIOS
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • NEMO/branches/2019/dev_r11351_fldread_with_XIOS

    • Property svn:externals

      •  

        old new  
        33^/utils/build/mk@HEAD         mk 
        44^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         5^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
         8 
         9# SETTE 
         10^/utils/CI/sette@13382        sette 

  • NEMO/branches/2019/dev_r11351_fldread_with_XIOS/src/OCE/SBC/sbcdcy.F90

    r10425 r13463  
    77   !!   NEMO    2.0  !  2006-02  (S. Masson, G. Madec)  adaptation to NEMO 
    88   !!           3.1  !  2009-07  (J.M. Molines)  adaptation to v3.1 
     9   !!           4.*  !  2019-10  (L. Brodeau)  nothing really new, but the routine 
     10   !!                ! "sbc_dcy_param" has been extracted from old function "sbc_dcy" 
     11   !!                ! => this allows the warm-layer param of COARE3* to know the time 
     12   !!                ! of dawn and dusk even if "ln_dm2dc=.false." (rdawn_dcy & rdusk_dcy 
     13   !!                ! are now public) 
    914   !!---------------------------------------------------------------------- 
    1015 
    … …  
    2227   IMPLICIT NONE 
    2328   PRIVATE 
    24     
     29 
    2530   INTEGER, PUBLIC ::   nday_qsr   !: day when parameters were computed 
    26     
     31 
    2732   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   raa , rbb  , rcc  , rab     ! diurnal cycle parameters 
    28    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   rtmd, rdawn, rdusk, rscal   !    -      -       - 
    29    
     33   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   rtmd, rscal   !    -      -       - 
     34   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: rdawn_dcy, rdusk_dcy   !    -      -       - 
     35 
    3036   PUBLIC   sbc_dcy        ! routine called by sbc 
    31  
     37   PUBLIC   sbc_dcy_param  ! routine used here and called by warm-layer parameterization (sbcblk_skin_coare*) 
     38 
     39   !! * Substitutions 
     40#  include "do_loop_substitute.h90" 
    3241   !!---------------------------------------------------------------------- 
    3342   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
    34    !! $Id$  
     43   !! $Id$ 
    3544   !! Software governed by the CeCILL license (see ./LICENSE) 
    3645   !!---------------------------------------------------------------------- 
    3746CONTAINS 
    3847 
    39       INTEGER FUNCTION sbc_dcy_alloc() 
    40          !!---------------------------------------------------------------------- 
    41          !!                ***  FUNCTION sbc_dcy_alloc  *** 
    42          !!---------------------------------------------------------------------- 
    43          ALLOCATE( raa (jpi,jpj) , rbb  (jpi,jpj) , rcc  (jpi,jpj) , rab  (jpi,jpj) ,     & 
    44             &      rtmd(jpi,jpj) , rdawn(jpi,jpj) , rdusk(jpi,jpj) , rscal(jpi,jpj) , STAT=sbc_dcy_alloc ) 
    45             ! 
    46          CALL mpp_sum ( 'sbcdcy', sbc_dcy_alloc ) 
    47          IF( sbc_dcy_alloc /= 0 )   CALL ctl_stop( 'STOP', 'sbc_dcy_alloc: failed to allocate arrays' ) 
    48       END FUNCTION sbc_dcy_alloc 
     48   INTEGER FUNCTION sbc_dcy_alloc() 
     49      !!---------------------------------------------------------------------- 
     50      !!                ***  FUNCTION sbc_dcy_alloc  *** 
     51      !!---------------------------------------------------------------------- 
     52      ALLOCATE( raa (jpi,jpj) , rbb  (jpi,jpj) , rcc  (jpi,jpj) , rab  (jpi,jpj) ,     & 
     53         &      rtmd(jpi,jpj) , rdawn_dcy(jpi,jpj) , rdusk_dcy(jpi,jpj) , rscal(jpi,jpj) , STAT=sbc_dcy_alloc ) 
     54      ! 
     55      CALL mpp_sum ( 'sbcdcy', sbc_dcy_alloc ) 
     56      IF( sbc_dcy_alloc /= 0 )   CALL ctl_stop( 'STOP', 'sbc_dcy_alloc: failed to allocate arrays' ) 
     57   END FUNCTION sbc_dcy_alloc 
    4958 
    5059 
    … …  
    6069      !! 
    6170      !! reference  : Bernie, DJ, E Guilyardi, G Madec, JM Slingo, and SJ Woolnough, 2007 
    62       !!              Impact of resolving the diurnal cycle in an ocean--atmosphere GCM.  
     71      !!              Impact of resolving the diurnal cycle in an ocean--atmosphere GCM. 
    6372      !!              Part 1: a diurnally forced OGCM. Climate Dynamics 29:6, 575-590. 
    6473      !!---------------------------------------------------------------------- 
    6574      LOGICAL , OPTIONAL          , INTENT(in) ::   l_mask    ! use the routine for night mask computation 
    66       REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pqsrin    ! input daily QSR flux  
     75      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pqsrin    ! input daily QSR flux 
    6776      REAL(wp), DIMENSION(jpi,jpj)             ::   zqsrout   ! output QSR flux with diurnal cycle 
    6877      !! 
    6978      INTEGER  ::   ji, jj                                       ! dummy loop indices 
    7079      INTEGER, DIMENSION(jpi,jpj) :: imask_night ! night mask 
    71       REAL(wp) ::   ztwopi, zinvtwopi, zconvrad  
    7280      REAL(wp) ::   zlo, zup, zlousd, zupusd 
    73       REAL(wp) ::   zdsws, zdecrad, ztx, zsin, zcos 
    74       REAL(wp) ::   ztmp, ztmp1, ztmp2, ztest 
     81      REAL(wp) ::   ztmp, ztmp1, ztmp2 
    7582      REAL(wp) ::   ztmpm, ztmpm1, ztmpm2 
    76       !---------------------------statement functions------------------------ 
    77       REAL(wp) ::   fintegral, pt1, pt2, paaa, pbbb, pccc        ! dummy statement function arguments 
    78       fintegral( pt1, pt2, paaa, pbbb, pccc ) =                         & 
    79          &   paaa * pt2 + zinvtwopi * pbbb * SIN(pccc + ztwopi * pt2)   & 
    80          & - paaa * pt1 - zinvtwopi * pbbb * SIN(pccc + ztwopi * pt1) 
    8183      !!--------------------------------------------------------------------- 
    8284      ! 
    8385      ! Initialization 
    8486      ! -------------- 
    85       ztwopi    = 2._wp * rpi 
    86       zinvtwopi = 1._wp / ztwopi 
    87       zconvrad  = ztwopi / 360._wp 
    88  
    8987      ! When are we during the day (from 0 to 1) 
    90       zlo = ( REAL(nsec_day, wp) - 0.5_wp * rdt ) / rday 
    91       zup = zlo + ( REAL(nn_fsbc, wp)     * rdt ) / rday 
    92       !                                           
    93       IF( nday_qsr == -1 ) THEN       ! first time step only   
     88      zlo = ( REAL(nsec_day, wp) - 0.5_wp * rn_Dt ) / rday 
     89      zup = zlo + ( REAL(nn_fsbc, wp)     * rn_Dt ) / rday 
     90      ! 
     91      IF( nday_qsr == -1 ) THEN       ! first time step only 
    9492         IF(lwp) THEN 
    9593            WRITE(numout,*) 
    … …  
    9896            WRITE(numout,*) 
    9997         ENDIF 
     98      ENDIF 
     99 
     100      ! Setting parameters for each new day: 
     101      CALL sbc_dcy_param() 
     102 
     103      !CALL iom_put( "rdusk_dcy", rdusk_dcy(:,:)*tmask(:,:,1) ) !LB 
     104      !CALL iom_put( "rdawn_dcy", rdawn_dcy(:,:)*tmask(:,:,1) ) !LB 
     105      !CALL iom_put( "rscal_dcy", rscal(:,:)*tmask(:,:,1) ) !LB 
     106 
     107 
     108      !     3. update qsr with the diurnal cycle 
     109      !     ------------------------------------ 
     110 
     111      imask_night(:,:) = 0 
     112      DO_2D( 1, 1, 1, 1 ) 
     113         ztmpm = 0._wp 
     114         IF( ABS(rab(ji,jj)) < 1. ) THEN         ! day duration is less than 24h 
     115            ! 
     116            IF( rdawn_dcy(ji,jj) < rdusk_dcy(ji,jj) ) THEN       ! day time in one part 
     117               zlousd = MAX(zlo, rdawn_dcy(ji,jj)) 
     118               zlousd = MIN(zlousd, zup) 
     119               zupusd = MIN(zup, rdusk_dcy(ji,jj)) 
     120               zupusd = MAX(zupusd, zlo) 
     121               ztmp = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
     122               zqsrout(ji,jj) = pqsrin(ji,jj) * ztmp * rscal(ji,jj) 
     123               ztmpm = zupusd - zlousd 
     124               IF( ztmpm .EQ. 0 ) imask_night(ji,jj) = 1 
     125               ! 
     126            ELSE                                         ! day time in two parts 
     127               zlousd = MIN(zlo, rdusk_dcy(ji,jj)) 
     128               zupusd = MIN(zup, rdusk_dcy(ji,jj)) 
     129               ztmp1 = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
     130               ztmpm1=zupusd-zlousd 
     131               zlousd = MAX(zlo, rdawn_dcy(ji,jj)) 
     132               zupusd = MAX(zup, rdawn_dcy(ji,jj)) 
     133               ztmp2 = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
     134               ztmpm2 =zupusd-zlousd 
     135               ztmp = ztmp1 + ztmp2 
     136               ztmpm = ztmpm1 + ztmpm2 
     137               zqsrout(ji,jj) = pqsrin(ji,jj) * ztmp * rscal(ji,jj) 
     138               IF(ztmpm .EQ. 0.) imask_night(ji,jj) = 1 
     139            ENDIF 
     140         ELSE                                   ! 24h light or 24h night 
     141            ! 
     142            IF( raa(ji,jj) > rbb(ji,jj) ) THEN           ! 24h day 
     143               ztmp = fintegral(zlo, zup, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
     144               zqsrout(ji,jj) = pqsrin(ji,jj) * ztmp * rscal(ji,jj) 
     145               imask_night(ji,jj) = 0 
     146               ! 
     147            ELSE                                         ! No day 
     148               zqsrout(ji,jj) = 0.0_wp 
     149               imask_night(ji,jj) = 1 
     150            ENDIF 
     151         ENDIF 
     152      END_2D 
     153      ! 
     154      IF( PRESENT(l_mask) .AND. l_mask ) THEN 
     155         zqsrout(:,:) = float(imask_night(:,:)) 
     156      ENDIF 
     157      ! 
     158   END FUNCTION sbc_dcy 
     159 
     160 
     161   SUBROUTINE sbc_dcy_param( ) 
     162      !! 
     163      INTEGER  ::   ji, jj                                       ! dummy loop indices 
     164      !INTEGER, DIMENSION(jpi,jpj) :: imask_night ! night mask 
     165      REAL(wp) ::   zdsws, zdecrad, ztx, zsin, zcos 
     166      REAL(wp) ::   ztmp, ztest 
     167      !---------------------------statement functions------------------------ 
     168      ! 
     169      IF( nday_qsr == -1 ) THEN       ! first time step only 
    100170         ! allocate sbcdcy arrays 
    101171         IF( sbc_dcy_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_dcy_alloc : unable to allocate arrays' ) 
    102172         ! Compute rcc needed to compute the time integral of the diurnal cycle 
    103          rcc(:,:) = zconvrad * glamt(:,:) - rpi 
     173         rcc(:,:) = rad * glamt(:,:) - rpi 
    104174         ! time of midday 
    105175         rtmd(:,:) = 0.5_wp - glamt(:,:) / 360._wp 
    … …  
    107177      ENDIF 
    108178 
    109       ! If this is a new day, we have to update the dawn, dusk and scaling function   
     179      ! If this is a new day, we have to update the dawn, dusk and scaling function 
    110180      !---------------------- 
    111      
    112       !     2.1 dawn and dusk   
    113  
    114       ! nday is the number of days since the beginning of the current month  
    115       IF( nday_qsr /= nday ) THEN  
     181 
     182      !     2.1 dawn and dusk 
     183 
     184      ! nday is the number of days since the beginning of the current month 
     185      IF( nday_qsr /= nday ) THEN 
    116186         ! save the day of the year and the daily mean of qsr 
    117          nday_qsr = nday  
    118          ! number of days since the previous winter solstice (supposed to be always 21 December)          
     187         nday_qsr = nday 
     188         ! number of days since the previous winter solstice (supposed to be always 21 December) 
    119189         zdsws = REAL(11 + nday_year, wp) 
    120190         ! declination of the earths orbit 
    121          zdecrad = (-23.5_wp * zconvrad) * COS( zdsws * ztwopi / REAL(nyear_len(1),wp) ) 
     191         zdecrad = (-23.5_wp * rad) * COS( zdsws * 2._wp*rpi / REAL(nyear_len(1),wp) ) 
    122192         ! Compute A and B needed to compute the time integral of the diurnal cycle 
    123193 
    124194         zsin = SIN( zdecrad )   ;   zcos = COS( zdecrad ) 
    125          DO jj = 1, jpj 
    126             DO ji = 1, jpi 
    127                ztmp = zconvrad * gphit(ji,jj) 
    128                raa(ji,jj) = SIN( ztmp ) * zsin 
    129                rbb(ji,jj) = COS( ztmp ) * zcos 
    130             END DO   
    131          END DO   
     195         DO_2D( 1, 1, 1, 1 ) 
     196            ztmp = rad * gphit(ji,jj) 
     197            raa(ji,jj) = SIN( ztmp ) * zsin 
     198            rbb(ji,jj) = COS( ztmp ) * zcos 
     199         END_2D 
    132200         ! Compute the time of dawn and dusk 
    133201 
    134          ! rab to test if the day time is equal to 0, less than 24h of full day         
     202         ! rab to test if the day time is equal to 0, less than 24h of full day 
    135203         rab(:,:) = -raa(:,:) / rbb(:,:) 
    136          DO jj = 1, jpj 
    137             DO ji = 1, jpi 
    138                IF ( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h 
    139          ! When is it night? 
    140                   ztx = zinvtwopi * (ACOS(rab(ji,jj)) - rcc(ji,jj)) 
    141                   ztest = -rbb(ji,jj) * SIN( rcc(ji,jj) + ztwopi * ztx ) 
    142          ! is it dawn or dusk? 
    143                   IF ( ztest > 0._wp ) THEN 
    144                      rdawn(ji,jj) = ztx 
    145                      rdusk(ji,jj) = rtmd(ji,jj) + ( rtmd(ji,jj) - rdawn(ji,jj) ) 
    146                   ELSE 
    147                      rdusk(ji,jj) = ztx 
    148                      rdawn(ji,jj) = rtmd(ji,jj) - ( rdusk(ji,jj) - rtmd(ji,jj) ) 
    149                   ENDIF 
     204         DO_2D( 1, 1, 1, 1 ) 
     205            IF( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h 
     206               ! When is it night? 
     207               ztx = 1._wp/(2._wp*rpi) * (ACOS(rab(ji,jj)) - rcc(ji,jj)) 
     208               ztest = -rbb(ji,jj) * SIN( rcc(ji,jj) + 2._wp*rpi * ztx ) 
     209               ! is it dawn or dusk? 
     210               IF( ztest > 0._wp ) THEN 
     211                  rdawn_dcy(ji,jj) = ztx 
     212                  rdusk_dcy(ji,jj) = rtmd(ji,jj) + ( rtmd(ji,jj) - rdawn_dcy(ji,jj) ) 
    150213               ELSE 
    151                   rdawn(ji,jj) = rtmd(ji,jj) + 0.5_wp 
    152                   rdusk(ji,jj) = rdawn(ji,jj) 
     214                  rdusk_dcy(ji,jj) = ztx 
     215                  rdawn_dcy(ji,jj) = rtmd(ji,jj) - ( rdusk_dcy(ji,jj) - rtmd(ji,jj) ) 
    153216               ENDIF 
    154              END DO   
    155          END DO   
    156          rdawn(:,:) = MOD( (rdawn(:,:) + 1._wp), 1._wp ) 
    157          rdusk(:,:) = MOD( (rdusk(:,:) + 1._wp), 1._wp ) 
     217            ELSE 
     218               rdawn_dcy(ji,jj) = rtmd(ji,jj) + 0.5_wp 
     219               rdusk_dcy(ji,jj) = rdawn_dcy(ji,jj) 
     220            ENDIF 
     221         END_2D 
     222         rdawn_dcy(:,:) = MOD( (rdawn_dcy(:,:) + 1._wp), 1._wp ) 
     223         rdusk_dcy(:,:) = MOD( (rdusk_dcy(:,:) + 1._wp), 1._wp ) 
    158224         !     2.2 Compute the scaling function: 
    159225         !         S* = the inverse of the time integral of the diurnal cycle from dawn to dusk 
    160226         !         Avoid possible infinite scaling factor, associated with very short daylight 
    161227         !         periods, by ignoring periods less than 1/1000th of a day (ticket #1040) 
    162          DO jj = 1, jpj 
    163             DO ji = 1, jpi 
    164                IF ( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h 
    165                   rscal(ji,jj) = 0.0_wp 
    166                   IF ( rdawn(ji,jj) < rdusk(ji,jj) ) THEN      ! day time in one part 
    167                      IF( (rdusk(ji,jj) - rdawn(ji,jj) ) .ge. 0.001_wp ) THEN 
    168                        rscal(ji,jj) = fintegral(rdawn(ji,jj), rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    169                        rscal(ji,jj) = 1._wp / rscal(ji,jj) 
    170                      ENDIF 
    171                   ELSE                                         ! day time in two parts 
    172                      IF( (rdusk(ji,jj) + (1._wp - rdawn(ji,jj)) ) .ge. 0.001_wp ) THEN 
    173                        rscal(ji,jj) = fintegral(0._wp, rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))   & 
    174                           &         + fintegral(rdawn(ji,jj), 1._wp, raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    175                        rscal(ji,jj) = 1. / rscal(ji,jj) 
    176                      ENDIF 
     228         DO_2D( 1, 1, 1, 1 ) 
     229            IF( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h 
     230               rscal(ji,jj) = 0.0_wp 
     231               IF( rdawn_dcy(ji,jj) < rdusk_dcy(ji,jj) ) THEN      ! day time in one part 
     232                  IF( (rdusk_dcy(ji,jj) - rdawn_dcy(ji,jj) ) .ge. 0.001_wp ) THEN 
     233                     rscal(ji,jj) = fintegral(rdawn_dcy(ji,jj), rdusk_dcy(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
     234                     rscal(ji,jj) = 1._wp / rscal(ji,jj) 
    177235                  ENDIF 
    178                ELSE 
    179                   IF ( raa(ji,jj) > rbb(ji,jj) ) THEN         ! 24h day 
    180                      rscal(ji,jj) = fintegral(0._wp, 1._wp, raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    181                      rscal(ji,jj) = 1._wp / rscal(ji,jj) 
    182                   ELSE                                          ! No day 
    183                      rscal(ji,jj) = 0.0_wp 
     236               ELSE                                         ! day time in two parts 
     237                  IF( (rdusk_dcy(ji,jj) + (1._wp - rdawn_dcy(ji,jj)) ) .ge. 0.001_wp ) THEN 
     238                     rscal(ji,jj) = fintegral(0._wp, rdusk_dcy(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))   & 
     239                        &         + fintegral(rdawn_dcy(ji,jj), 1._wp, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
     240                     rscal(ji,jj) = 1. / rscal(ji,jj) 
    184241                  ENDIF 
    185242               ENDIF 
    186             END DO   
    187          END DO   
     243            ELSE 
     244               IF( raa(ji,jj) > rbb(ji,jj) ) THEN         ! 24h day 
     245                  rscal(ji,jj) = fintegral(0._wp, 1._wp, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
     246                  rscal(ji,jj) = 1._wp / rscal(ji,jj) 
     247               ELSE                                          ! No day 
     248                  rscal(ji,jj) = 0.0_wp 
     249               ENDIF 
     250            ENDIF 
     251         END_2D 
    188252         ! 
    189          ztmp = rday / ( rdt * REAL(nn_fsbc, wp) ) 
     253         ztmp = rday / ( rn_Dt * REAL(nn_fsbc, wp) ) 
    190254         rscal(:,:) = rscal(:,:) * ztmp 
    191255         ! 
    192       ENDIF  
    193          !     3. update qsr with the diurnal cycle 
    194          !     ------------------------------------ 
    195  
    196       imask_night(:,:) = 0 
    197       DO jj = 1, jpj 
    198          DO ji = 1, jpi 
    199             ztmpm = 0._wp 
    200             IF( ABS(rab(ji,jj)) < 1. ) THEN         ! day duration is less than 24h 
    201                ! 
    202                IF( rdawn(ji,jj) < rdusk(ji,jj) ) THEN       ! day time in one part 
    203                   zlousd = MAX(zlo, rdawn(ji,jj)) 
    204                   zlousd = MIN(zlousd, zup) 
    205                   zupusd = MIN(zup, rdusk(ji,jj)) 
    206                   zupusd = MAX(zupusd, zlo) 
    207                   ztmp = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    208                   zqsrout(ji,jj) = pqsrin(ji,jj) * ztmp * rscal(ji,jj) 
    209                   ztmpm = zupusd - zlousd 
    210                   IF ( ztmpm .EQ. 0 ) imask_night(ji,jj) = 1 
    211                   ! 
    212                ELSE                                         ! day time in two parts 
    213                   zlousd = MIN(zlo, rdusk(ji,jj)) 
    214                   zupusd = MIN(zup, rdusk(ji,jj)) 
    215                   ztmp1 = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    216                   ztmpm1=zupusd-zlousd 
    217                   zlousd = MAX(zlo, rdawn(ji,jj)) 
    218                   zupusd = MAX(zup, rdawn(ji,jj)) 
    219                   ztmp2 = fintegral(zlousd, zupusd, raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    220                   ztmpm2 =zupusd-zlousd 
    221                   ztmp = ztmp1 + ztmp2 
    222                   ztmpm = ztmpm1 + ztmpm2 
    223                   zqsrout(ji,jj) = pqsrin(ji,jj) * ztmp * rscal(ji,jj) 
    224                   IF (ztmpm .EQ. 0.) imask_night(ji,jj) = 1 
    225                ENDIF 
    226             ELSE                                   ! 24h light or 24h night 
    227                ! 
    228                IF( raa(ji,jj) > rbb(ji,jj) ) THEN           ! 24h day 
    229                   ztmp = fintegral(zlo, zup, raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    230                   zqsrout(ji,jj) = pqsrin(ji,jj) * ztmp * rscal(ji,jj) 
    231                   imask_night(ji,jj) = 0 
    232                   ! 
    233                ELSE                                         ! No day 
    234                   zqsrout(ji,jj) = 0.0_wp 
    235                   imask_night(ji,jj) = 1 
    236                ENDIF 
    237             ENDIF 
    238          END DO   
    239       END DO   
    240       ! 
    241       IF( PRESENT(l_mask) .AND. l_mask ) THEN 
    242          zqsrout(:,:) = float(imask_night(:,:)) 
    243       ENDIF 
    244       ! 
    245    END FUNCTION sbc_dcy 
     256      ENDIF !IF( nday_qsr /= nday ) 
     257      ! 
     258   END SUBROUTINE sbc_dcy_param 
     259 
     260 
     261   FUNCTION fintegral( pt1, pt2, paaa, pbbb, pccc ) 
     262      REAL(wp), INTENT(in) :: pt1, pt2, paaa, pbbb, pccc 
     263      REAL(wp) :: fintegral 
     264      fintegral =   paaa * pt2 + 1._wp/(2._wp*rpi) * pbbb * SIN(pccc + 2._wp*rpi*pt2)   & 
     265         &        - paaa * pt1 - 1._wp/(2._wp*rpi) * pbbb * SIN(pccc + 2._wp*rpi*pt1) 
     266   END FUNCTION fintegral 
    246267 
    247268   !!====================================================================== 
Note: See TracChangeset for help on using the changeset viewer.