source: NEMO/releases/r4.0/r4.0-HEAD/src/OCE/TRA/traqsr.F90 @ 14715

Last change on this file since 14715 was 14715, checked in by clem, 6 months ago

solve ticket #2641

  • Property svn:keywords set to Id
File size: 22.5 KB
Line 
1MODULE traqsr
2   !!======================================================================
3   !!                       ***  MODULE  traqsr  ***
4   !! Ocean physics:   solar radiation penetration in the top ocean levels
5   !!======================================================================
6   !! History :  OPA  !  1990-10  (B. Blanke)  Original code
7   !!            7.0  !  1991-11  (G. Madec)
8   !!                 !  1996-01  (G. Madec)  s-coordinates
9   !!   NEMO     1.0  !  2002-06  (G. Madec)  F90: Free form and module
10   !!             -   !  2005-11  (G. Madec) zco, zps, sco coordinate
11   !!            3.2  !  2009-04  (G. Madec & NEMO team)
12   !!            3.6  !  2012-05  (C. Rousset) store attenuation coef for use in ice model
13   !!            3.6  !  2015-12  (O. Aumont, J. Jouanno, C. Ethe) use vertical profile of chlorophyll
14   !!            3.7  !  2015-11  (G. Madec, A. Coward)  remove optimisation for fix volume
15   !!----------------------------------------------------------------------
16
17   !!----------------------------------------------------------------------
18   !!   tra_qsr       : temperature trend due to the penetration of solar radiation
19   !!   tra_qsr_init  : initialization of the qsr penetration
20   !!----------------------------------------------------------------------
21   USE oce            ! ocean dynamics and active tracers
22   USE phycst         ! physical constants
23   USE dom_oce        ! ocean space and time domain
24   USE sbc_oce        ! surface boundary condition: ocean
25   USE trc_oce        ! share SMS/Ocean variables
26   USE trd_oce        ! trends: ocean variables
27   USE trdtra         ! trends manager: tracers
28   !
29   USE in_out_manager ! I/O manager
30   USE prtctl         ! Print control
31   USE iom            ! I/O library
32   USE fldread        ! read input fields
33   USE restart        ! ocean restart
34   USE lib_mpp        ! MPP library
35   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
36   USE timing         ! Timing
37
38   IMPLICIT NONE
39   PRIVATE
40
41   PUBLIC   tra_qsr       ! routine called by step.F90 (ln_traqsr=T)
42   PUBLIC   tra_qsr_init  ! routine called by nemogcm.F90
43
44   !                                 !!* Namelist namtra_qsr: penetrative solar radiation
45   LOGICAL , PUBLIC ::   ln_traqsr    !: light absorption (qsr) flag
46   LOGICAL , PUBLIC ::   ln_qsr_rgb   !: Red-Green-Blue light absorption flag 
47   LOGICAL , PUBLIC ::   ln_qsr_2bd   !: 2 band         light absorption flag
48   LOGICAL , PUBLIC ::   ln_qsr_bio   !: bio-model      light absorption flag
49   INTEGER , PUBLIC ::   nn_chldta    !: use Chlorophyll data (=1) or not (=0)
50   REAL(wp), PUBLIC ::   rn_abs       !: fraction absorbed in the very near surface (RGB & 2 bands)
51   REAL(wp), PUBLIC ::   rn_si0       !: very near surface depth of extinction      (RGB & 2 bands)
52   REAL(wp), PUBLIC ::   rn_si1       !: deepest depth of extinction (water type I)       (2 bands)
53   !
54   INTEGER , PUBLIC ::   nksr         !: levels below which the light cannot penetrate (depth larger than 391 m)
55 
56   INTEGER, PARAMETER ::   np_RGB  = 1   ! R-G-B     light penetration with constant Chlorophyll
57   INTEGER, PARAMETER ::   np_RGBc = 2   ! R-G-B     light penetration with Chlorophyll data
58   INTEGER, PARAMETER ::   np_2BD  = 3   ! 2 bands   light penetration
59   INTEGER, PARAMETER ::   np_BIO  = 4   ! bio-model light penetration
60   !
61   INTEGER  ::   nqsr    ! user choice of the type of light penetration
62   REAL(wp) ::   xsi0r   ! inverse of rn_si0
63   REAL(wp) ::   xsi1r   ! inverse of rn_si1
64   !
65   REAL(wp) , PUBLIC, DIMENSION(3,61)   ::   rkrgb    ! tabulated attenuation coefficients for RGB absorption
66   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_chl   ! structure of input Chl (file informations, fields read)
67
68   !! * Substitutions
69#  include "vectopt_loop_substitute.h90"
70   !!----------------------------------------------------------------------
71   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
72   !! $Id$
73   !! Software governed by the CeCILL license (see ./LICENSE)
74   !!----------------------------------------------------------------------
75CONTAINS
76
77   SUBROUTINE tra_qsr( kt )
78      !!----------------------------------------------------------------------
79      !!                  ***  ROUTINE tra_qsr  ***
80      !!
81      !! ** Purpose :   Compute the temperature trend due to the solar radiation
82      !!              penetration and add it to the general temperature trend.
83      !!
84      !! ** Method  : The profile of the solar radiation within the ocean is defined
85      !!      through 2 wavebands (rn_si0,rn_si1) or 3 wavebands (RGB) and a ratio rn_abs
86      !!      Considering the 2 wavebands case:
87      !!         I(k) = Qsr*( rn_abs*EXP(z(k)/rn_si0) + (1.-rn_abs)*EXP(z(k)/rn_si1) )
88      !!         The temperature trend associated with the solar radiation penetration
89      !!         is given by : zta = 1/e3t dk[ I ] / (rau0*Cp)
90      !!         At the bottom, boudary condition for the radiation is no flux :
91      !!      all heat which has not been absorbed in the above levels is put
92      !!      in the last ocean level.
93      !!         The computation is only done down to the level where
94      !!      I(k) < 1.e-15 W/m2 (i.e. over the top nksr levels) .
95      !!
96      !! ** Action  : - update ta with the penetrative solar radiation trend
97      !!              - send  trend for further diagnostics (l_trdtra=T)
98      !!
99      !! Reference  : Jerlov, N. G., 1968 Optical Oceanography, Elsevier, 194pp.
100      !!              Lengaigne et al. 2007, Clim. Dyn., V28, 5, 503-516.
101      !!              Morel, A. et Berthon, JF, 1989, Limnol Oceanogr 34(8), 1545-1562
102      !!----------------------------------------------------------------------
103      INTEGER, INTENT(in) ::   kt     ! ocean time-step
104      !
105      INTEGER  ::   ji, jj, jk               ! dummy loop indices
106      INTEGER  ::   irgb                     ! local integers
107      REAL(wp) ::   zchl, zcoef, z1_2        ! local scalars
108      REAL(wp) ::   zc0 , zc1 , zc2 , zc3    !    -         -
109      REAL(wp) ::   zzc0, zzc1, zzc2, zzc3   !    -         -
110      REAL(wp) ::   zz0 , zz1                !    -         -
111      REAL(wp) ::   zCb, zCmax, zze, zpsi, zpsimax, zdelpsi, zCtot, zCze
112      REAL(wp) ::   zlogc, zlogc2, zlogc3 
113      REAL(wp), ALLOCATABLE, DIMENSION(:,:)   :: zekb, zekg, zekr
114      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ze0, ze1, ze2, ze3, zea, ztrdt
115      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zetot, zchl3d
116      !!----------------------------------------------------------------------
117      !
118      IF( ln_timing )   CALL timing_start('tra_qsr')
119      !
120      IF( kt == nit000 ) THEN
121         IF(lwp) WRITE(numout,*)
122         IF(lwp) WRITE(numout,*) 'tra_qsr : penetration of the surface solar radiation'
123         IF(lwp) WRITE(numout,*) '~~~~~~~'
124      ENDIF
125      !
126      IF( l_trdtra ) THEN      ! trends diagnostic: save the input temperature trend
127         ALLOCATE( ztrdt(jpi,jpj,jpk) ) 
128         ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
129      ENDIF
130      !
131      !                         !-----------------------------------!
132      !                         !  before qsr induced heat content  !
133      !                         !-----------------------------------!
134      IF( kt == nit000 ) THEN          !==  1st time step  ==!
135!!gm case neuler  not taken into account....
136         IF( ln_rstart .AND. iom_varid( numror, 'qsr_hc_b', ldstop = .FALSE. ) > 0 ) THEN    ! read in restart
137            IF(lwp) WRITE(numout,*) '          nit000-1 qsr tracer content forcing field read in the restart file'
138            z1_2 = 0.5_wp
139            CALL iom_get( numror, jpdom_autoglo, 'qsr_hc_b', qsr_hc_b, ldxios = lrxios )   ! before heat content trend due to Qsr flux
140         ELSE                                           ! No restart or restart not found: Euler forward time stepping
141            z1_2 = 1._wp
142            qsr_hc_b(:,:,:) = 0._wp
143         ENDIF
144      ELSE                             !==  Swap of qsr heat content  ==!
145         z1_2 = 0.5_wp
146         qsr_hc_b(:,:,:) = qsr_hc(:,:,:)
147      ENDIF
148      !
149      !                         !--------------------------------!
150      SELECT CASE( nqsr )       !  now qsr induced heat content  !
151      !                         !--------------------------------!
152      !
153      CASE( np_BIO )                   !==  bio-model fluxes  ==!
154         !
155         DO jk = 1, nksr
156            qsr_hc(:,:,jk) = r1_rau0_rcp * ( etot3(:,:,jk) - etot3(:,:,jk+1) )
157         END DO
158         !
159      CASE( np_RGB , np_RGBc )         !==  R-G-B fluxes  ==!
160         !
161         ALLOCATE( zekb(jpi,jpj)     , zekg(jpi,jpj)     , zekr  (jpi,jpj)     , &
162            &      ze0 (jpi,jpj,jpk) , ze1 (jpi,jpj,jpk) , ze2   (jpi,jpj,jpk) , &
163            &      ze3 (jpi,jpj,jpk) , zea (jpi,jpj,jpk) , zchl3d(jpi,jpj,jpk)   ) 
164         !
165         IF( nqsr == np_RGBc ) THEN          !*  Variable Chlorophyll
166            CALL fld_read( kt, 1, sf_chl )         ! Read Chl data and provides it at the current time step
167            DO jk = 1, nksr + 1
168               DO jj = 2, jpjm1                       ! Separation in R-G-B depending of the surface Chl
169                  DO ji = fs_2, fs_jpim1
170                     zchl    = MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) )
171                     zCtot   = 40.6  * zchl**0.459
172                     zze     = 568.2 * zCtot**(-0.746)
173                     IF( zze > 102. ) zze = 200.0 * zCtot**(-0.293)
174                     zpsi    = gdepw_n(ji,jj,jk) / zze
175                     !
176                     zlogc   = LOG( zchl )
177                     zlogc2  = zlogc * zlogc
178                     zlogc3  = zlogc * zlogc * zlogc
179                     zCb     = 0.768 + 0.087 * zlogc - 0.179 * zlogc2 - 0.025 * zlogc3
180                     zCmax   = 0.299 - 0.289 * zlogc + 0.579 * zlogc2
181                     zpsimax = 0.6   - 0.640 * zlogc + 0.021 * zlogc2 + 0.115 * zlogc3
182                     zdelpsi = 0.710 + 0.159 * zlogc + 0.021 * zlogc2
183                     zCze    = 1.12  * (zchl)**0.803 
184                     !
185                     zchl3d(ji,jj,jk) = zCze * ( zCb + zCmax * EXP( -( (zpsi - zpsimax) / zdelpsi )**2 ) )
186                  END DO
187                  !
188               END DO
189            END DO
190         ELSE                                !* constant chrlorophyll
191           DO jk = 1, nksr + 1
192              zchl3d(:,:,jk) = 0.05 
193            ENDDO
194         ENDIF
195         !
196         zcoef  = ( 1. - rn_abs ) / 3._wp    !* surface equi-partition in R-G-B
197         DO jj = 2, jpjm1
198            DO ji = fs_2, fs_jpim1
199               ze0(ji,jj,1) = rn_abs * qsr(ji,jj)
200               ze1(ji,jj,1) = zcoef  * qsr(ji,jj)
201               ze2(ji,jj,1) = zcoef  * qsr(ji,jj)
202               ze3(ji,jj,1) = zcoef  * qsr(ji,jj)
203               zea(ji,jj,1) =          qsr(ji,jj)
204            END DO
205         END DO
206         !
207         DO jk = 2, nksr+1                   !* interior equi-partition in R-G-B depending of vertical profile of Chl
208            DO jj = 2, jpjm1
209               DO ji = fs_2, fs_jpim1
210                  zchl = MIN( 10. , MAX( 0.03, zchl3d(ji,jj,jk) ) )
211                  irgb = NINT( 41 + 20.*LOG10(zchl) + 1.e-15 )
212                  zekb(ji,jj) = rkrgb(1,irgb)
213                  zekg(ji,jj) = rkrgb(2,irgb)
214                  zekr(ji,jj) = rkrgb(3,irgb)
215               END DO
216            END DO
217
218            DO jj = 2, jpjm1
219               DO ji = fs_2, fs_jpim1
220                  zc0 = ze0(ji,jj,jk-1) * EXP( - e3t_n(ji,jj,jk-1) * xsi0r       )
221                  zc1 = ze1(ji,jj,jk-1) * EXP( - e3t_n(ji,jj,jk-1) * zekb(ji,jj) )
222                  zc2 = ze2(ji,jj,jk-1) * EXP( - e3t_n(ji,jj,jk-1) * zekg(ji,jj) )
223                  zc3 = ze3(ji,jj,jk-1) * EXP( - e3t_n(ji,jj,jk-1) * zekr(ji,jj) )
224                  ze0(ji,jj,jk) = zc0
225                  ze1(ji,jj,jk) = zc1
226                  ze2(ji,jj,jk) = zc2
227                  ze3(ji,jj,jk) = zc3
228                  zea(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * wmask(ji,jj,jk)
229               END DO
230            END DO
231         END DO
232         !
233         DO jk = 1, nksr                     !* now qsr induced heat content
234            DO jj = 2, jpjm1
235               DO ji = fs_2, fs_jpim1
236                  qsr_hc(ji,jj,jk) = r1_rau0_rcp * ( zea(ji,jj,jk) - zea(ji,jj,jk+1) )
237               END DO
238            END DO
239         END DO
240         !
241         DEALLOCATE( zekb , zekg , zekr , ze0 , ze1 , ze2 , ze3 , zea , zchl3d ) 
242         !
243      CASE( np_2BD  )            !==  2-bands fluxes  ==!
244         !
245         zz0 =        rn_abs   * r1_rau0_rcp      ! surface equi-partition in 2-bands
246         zz1 = ( 1. - rn_abs ) * r1_rau0_rcp
247         DO jk = 1, nksr                          ! solar heat absorbed at T-point in the top 400m
248            DO jj = 2, jpjm1
249               DO ji = fs_2, fs_jpim1
250                  zc0 = zz0 * EXP( -gdepw_n(ji,jj,jk  )*xsi0r ) + zz1 * EXP( -gdepw_n(ji,jj,jk  )*xsi1r )
251                  zc1 = zz0 * EXP( -gdepw_n(ji,jj,jk+1)*xsi0r ) + zz1 * EXP( -gdepw_n(ji,jj,jk+1)*xsi1r )
252                  qsr_hc(ji,jj,jk) = qsr(ji,jj) * ( zc0 * wmask(ji,jj,jk) - zc1 * wmask(ji,jj,jk+1) ) 
253               END DO
254            END DO
255         END DO
256         !
257      END SELECT
258      !
259      !                          !-----------------------------!
260      DO jk = 1, nksr            !  update to the temp. trend  !
261         DO jj = 2, jpjm1        !-----------------------------!
262            DO ji = fs_2, fs_jpim1   ! vector opt.
263               tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)   &
264                  &                 + z1_2 * ( qsr_hc_b(ji,jj,jk) + qsr_hc(ji,jj,jk) ) / e3t_n(ji,jj,jk)
265            END DO
266         END DO
267      END DO
268      !
269      ! sea-ice: store the 1st ocean level attenuation coefficient
270      DO jj = 2, jpjm1 
271         DO ji = fs_2, fs_jpim1   ! vector opt.
272            IF( qsr(ji,jj) /= 0._wp ) THEN   ;   fraqsr_1lev(ji,jj) = qsr_hc(ji,jj,1) / ( r1_rau0_rcp * qsr(ji,jj) )
273            ELSE                             ;   fraqsr_1lev(ji,jj) = 1._wp
274            ENDIF
275         END DO
276      END DO
277      CALL lbc_lnk( 'traqsr', fraqsr_1lev(:,:), 'T', 1._wp )
278      !
279      IF( iom_use('qsr3d') ) THEN      ! output the shortwave Radiation distribution
280         ALLOCATE( zetot(jpi,jpj,jpk) )
281         zetot(:,:,nksr+1:jpk) = 0._wp     ! below ~400m set to zero
282         DO jk = nksr, 1, -1
283            DO jj = 2, jpjm1
284               DO ji = fs_2, fs_jpim1
285                  zetot(ji,jj,jk) = zetot(ji,jj,jk+1) + qsr_hc(ji,jj,jk) * rau0_rcp
286               ENDDO
287            ENDDO
288         END DO
289         CALL lbc_lnk( 'traqsr', zetot, 'T', 1._wp )         
290         CALL iom_put( 'qsr3d', zetot )   ! 3D distribution of shortwave Radiation
291         DEALLOCATE( zetot ) 
292      ENDIF
293      !
294      IF( lrst_oce ) THEN     ! write in the ocean restart file
295         IF( lwxios ) CALL iom_swap(      cwxios_context          )
296         CALL iom_rstput( kt, nitrst, numrow, 'qsr_hc_b'   , qsr_hc     , ldxios = lwxios )
297         CALL iom_rstput( kt, nitrst, numrow, 'fraqsr_1lev', fraqsr_1lev, ldxios = lwxios ) 
298         IF( lwxios ) CALL iom_swap(      cxios_context          )
299      ENDIF
300      !
301      IF( l_trdtra ) THEN     ! qsr tracers trends saved for diagnostics
302         ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
303         CALL trd_tra( kt, 'TRA', jp_tem, jptra_qsr, ztrdt )
304         DEALLOCATE( ztrdt ) 
305      ENDIF
306      !                       ! print mean trends (used for debugging)
307      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' qsr  - Ta: ', mask1=tmask, clinfo3='tra-ta' )
308      !
309      IF( ln_timing )   CALL timing_stop('tra_qsr')
310      !
311   END SUBROUTINE tra_qsr
312
313
314   SUBROUTINE tra_qsr_init
315      !!----------------------------------------------------------------------
316      !!                  ***  ROUTINE tra_qsr_init  ***
317      !!
318      !! ** Purpose :   Initialization for the penetrative solar radiation
319      !!
320      !! ** Method  :   The profile of solar radiation within the ocean is set
321      !!      from two length scale of penetration (rn_si0,rn_si1) and a ratio
322      !!      (rn_abs). These parameters are read in the namtra_qsr namelist. The
323      !!      default values correspond to clear water (type I in Jerlov'
324      !!      (1968) classification.
325      !!         called by tra_qsr at the first timestep (nit000)
326      !!
327      !! ** Action  : - initialize rn_si0, rn_si1 and rn_abs
328      !!
329      !! Reference : Jerlov, N. G., 1968 Optical Oceanography, Elsevier, 194pp.
330      !!----------------------------------------------------------------------
331      INTEGER  ::   ji, jj, jk                  ! dummy loop indices
332      INTEGER  ::   ios, irgb, ierror, ioptio   ! local integer
333      REAL(wp) ::   zz0, zc0 , zc1, zcoef      ! local scalars
334      REAL(wp) ::   zz1, zc2 , zc3, zchl       !   -      -
335      !
336      CHARACTER(len=100) ::   cn_dir   ! Root directory for location of ssr files
337      TYPE(FLD_N)        ::   sn_chl   ! informations about the chlorofyl field to be read
338      !!
339      NAMELIST/namtra_qsr/  sn_chl, cn_dir, ln_qsr_rgb, ln_qsr_2bd, ln_qsr_bio,  &
340         &                  nn_chldta, rn_abs, rn_si0, rn_si1
341      !!----------------------------------------------------------------------
342      !
343      REWIND( numnam_ref )              ! Namelist namtra_qsr in reference     namelist
344      READ  ( numnam_ref, namtra_qsr, IOSTAT = ios, ERR = 901)
345901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtra_qsr in reference namelist' )
346      !
347      REWIND( numnam_cfg )              ! Namelist namtra_qsr in configuration namelist
348      READ  ( numnam_cfg, namtra_qsr, IOSTAT = ios, ERR = 902 )
349902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namtra_qsr in configuration namelist' )
350      IF(lwm) WRITE ( numond, namtra_qsr )
351      !
352      IF(lwp) THEN                ! control print
353         WRITE(numout,*)
354         WRITE(numout,*) 'tra_qsr_init : penetration of the surface solar radiation'
355         WRITE(numout,*) '~~~~~~~~~~~~'
356         WRITE(numout,*) '   Namelist namtra_qsr : set the parameter of penetration'
357         WRITE(numout,*) '      RGB (Red-Green-Blue) light penetration       ln_qsr_rgb = ', ln_qsr_rgb
358         WRITE(numout,*) '      2 band               light penetration       ln_qsr_2bd = ', ln_qsr_2bd
359         WRITE(numout,*) '      bio-model            light penetration       ln_qsr_bio = ', ln_qsr_bio
360         WRITE(numout,*) '      RGB : Chl data (=1) or cst value (=0)        nn_chldta  = ', nn_chldta
361         WRITE(numout,*) '      RGB & 2 bands: fraction of light (rn_si1)    rn_abs     = ', rn_abs
362         WRITE(numout,*) '      RGB & 2 bands: shortess depth of extinction  rn_si0     = ', rn_si0
363         WRITE(numout,*) '      2 bands: longest depth of extinction         rn_si1     = ', rn_si1
364         WRITE(numout,*)
365      ENDIF
366      !
367      ioptio = 0                    ! Parameter control
368      IF( ln_qsr_rgb  )   ioptio = ioptio + 1
369      IF( ln_qsr_2bd  )   ioptio = ioptio + 1
370      IF( ln_qsr_bio  )   ioptio = ioptio + 1
371      !
372      IF( ioptio /= 1 )   CALL ctl_stop( 'Choose ONE type of light penetration in namelist namtra_qsr',  &
373         &                               ' 2 bands, 3 RGB bands or bio-model light penetration' )
374      !
375      IF( ln_qsr_rgb .AND. nn_chldta == 0 )   nqsr = np_RGB 
376      IF( ln_qsr_rgb .AND. nn_chldta == 1 )   nqsr = np_RGBc
377      IF( ln_qsr_2bd                      )   nqsr = np_2BD
378      IF( ln_qsr_bio                      )   nqsr = np_BIO
379      !
380      !                             ! Initialisation
381      xsi0r = 1._wp / rn_si0
382      xsi1r = 1._wp / rn_si1
383      !
384      SELECT CASE( nqsr )
385      !                               
386      CASE( np_RGB , np_RGBc )         !==  Red-Green-Blue light penetration  ==!
387         !                             
388         IF(lwp)   WRITE(numout,*) '   ==>>>   R-G-B   light penetration '
389         !
390         CALL trc_oce_rgb( rkrgb )                 ! tabulated attenuation coef.
391         !                                   
392         nksr = trc_oce_ext_lev( r_si2, 33._wp )   ! level of light extinction
393         !
394         IF(lwp) WRITE(numout,*) '        level of light extinction = ', nksr, ' ref depth = ', gdepw_1d(nksr+1), ' m'
395         !
396         IF( nqsr == np_RGBc ) THEN                ! Chl data : set sf_chl structure
397            IF(lwp) WRITE(numout,*) '   ==>>>   Chlorophyll read in a file'
398            ALLOCATE( sf_chl(1), STAT=ierror )
399            IF( ierror > 0 ) THEN
400               CALL ctl_stop( 'tra_qsr_init: unable to allocate sf_chl structure' )   ;   RETURN
401            ENDIF
402            ALLOCATE( sf_chl(1)%fnow(jpi,jpj,1)   )
403            IF( sn_chl%ln_tint )   ALLOCATE( sf_chl(1)%fdta(jpi,jpj,1,2) )
404            !                                        ! fill sf_chl with sn_chl and control print
405            CALL fld_fill( sf_chl, (/ sn_chl /), cn_dir, 'tra_qsr_init',   &
406               &           'Solar penetration function of read chlorophyll', 'namtra_qsr' , no_print )
407         ENDIF
408         IF( nqsr == np_RGB ) THEN                 ! constant Chl
409            IF(lwp) WRITE(numout,*) '   ==>>>   Constant Chlorophyll concentration = 0.05'
410         ENDIF
411         !
412      CASE( np_2BD )                   !==  2 bands light penetration  ==!
413         !
414         IF(lwp)  WRITE(numout,*) '   ==>>>   2 bands light penetration'
415         !
416         nksr = trc_oce_ext_lev( rn_si1, 100._wp )    ! level of light extinction
417         IF(lwp) WRITE(numout,*) '        level of light extinction = ', nksr, ' ref depth = ', gdepw_1d(nksr+1), ' m'
418         !
419      CASE( np_BIO )                   !==  BIO light penetration  ==!
420         !
421         IF(lwp) WRITE(numout,*) '   ==>>>   bio-model light penetration'
422         IF( .NOT.lk_top )   CALL ctl_stop( 'No bio model : ln_qsr_bio = true impossible ' )
423         !                                   
424         CALL trc_oce_rgb( rkrgb )                 ! tabulated attenuation coef.
425         !                                   
426         nksr = trc_oce_ext_lev( r_si2, 33._wp )   ! level of light extinction
427         !
428         IF(lwp) WRITE(numout,*) '        level of light extinction = ', nksr, ' ref depth = ', gdepw_1d(nksr+1), ' m'
429         !
430      END SELECT
431      !
432      qsr_hc(:,:,:) = 0._wp     ! now qsr heat content set to zero where it will not be computed
433      !
434      ! 1st ocean level attenuation coefficient (used in sbcssm)
435      IF( iom_varid( numror, 'fraqsr_1lev', ldstop = .FALSE. ) > 0 ) THEN
436         CALL iom_get( numror, jpdom_autoglo, 'fraqsr_1lev'  , fraqsr_1lev, ldxios = lrxios  )
437      ELSE
438         fraqsr_1lev(:,:) = 1._wp   ! default : no penetration
439      ENDIF
440      !
441      IF( lwxios ) THEN
442         CALL iom_set_rstw_var_active('qsr_hc_b')
443         CALL iom_set_rstw_var_active('fraqsr_1lev')
444      ENDIF
445      !
446   END SUBROUTINE tra_qsr_init
447
448   !!======================================================================
449END MODULE traqsr
Note: See TracBrowser for help on using the repository browser.