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.
asmlogchlbal_ersem.F90 in branches/UKMO/dev_r5518_v3.6_asm_nemovar_community_logchlbal/NEMOGCM/NEMO/OPA_SRC/ASM – NEMO

source: branches/UKMO/dev_r5518_v3.6_asm_nemovar_community_logchlbal/NEMOGCM/NEMO/OPA_SRC/ASM/asmlogchlbal_ersem.F90 @ 8438

Last change on this file since 8438 was 8438, checked in by dford, 7 years ago

Initial implementation of Jozef S's code to do multivariate balancing with ERSEM.
File size: 17.6 KB
Line 
1MODULE asmlogchlbal_ersem
2   !!======================================================================
3   !!                       ***  MODULE asmlogchlbal_ersem  ***
4   !! Calculate increments to ERSEM based on surface logchl increments
5   !!======================================================================
6   !! History :  3.6  ! 2016-09  (D. Ford)     Original code
7   !!----------------------------------------------------------------------
8#if defined key_asminc && defined key_fabm
9   !!----------------------------------------------------------------------
10   !!   'key_asminc'         : assimilation increment interface
11   !!   'key_fabm'           : FABM-ERSEM coupling
12   !!----------------------------------------------------------------------
13   !!   asm_logchl_bal_ersem : routine to calculate increments to ERSEM
14   !!----------------------------------------------------------------------
15   USE par_kind,    ONLY: wp             ! kind parameters
16   USE par_oce,     ONLY: jpi, jpj, jpk  ! domain array sizes
17   USE dom_oce,     ONLY: gdepw_n        ! domain information
18   USE zdfmxl                            ! mixed layer depth
19   USE iom                               ! i/o
20   USE trc,         ONLY: trn            ! ERSEM state variables
21   USE par_fabm                          ! FABM parameters
22   USE par_trc,     ONLY: jptra          ! Tracer parameters
23
24
25   IMPLICIT NONE
26   PRIVATE                   
27
28   PUBLIC asm_logchl_bal_ersem
29
30   REAL(wp), PARAMETER :: tiny = 1e-20
31
32CONTAINS
33    SUBROUTINE asm_logchl_bal_ersem( ld_logchlpftinc, npfts, mld_choice_bgc, &
34      &                             k_maxchlinc, logchl_bkginc, logchl_balinc )                     
35
36                                                                       
37                                         
38 !! INTEGER, INTENT(in   ) :: npfts
39 !!REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj) :: inc_log_chlTot
40
41      LOGICAL,  INTENT(in   )                           :: ld_logchlpftinc
42      INTEGER,  INTENT(in   )                           :: npfts
43      INTEGER,  INTENT(in   )                           :: mld_choice_bgc
44      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,npfts) :: logchl_bkginc
45      REAL(wp), INTENT(in   )                               :: k_maxchlinc
46      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj,jpk,jptra)   :: logchl_balinc
47
48      !!
49 
50 INTEGER                      :: ji, jj, jk, jp
51 INTEGER                      :: jkmax
52 REAL(wp)                     :: chl_tot, chl_inc
53 REAL(wp), DIMENSION(jpi,jpj) :: zmld
54 REAL(wp), DIMENSION(jpi,jpj,npfts,jptra) :: cc_chl_field  !! CROSS-COV BETWEEN SPECIFIC SIZE CLASS AND GIVEN BIOGEOCHEM VAR
55 REAL(wp), DIMENSION(jptra) :: threshold
56 
57
58
59         ! STEP 1: TRANSFORM ASSIMILATED LOG-FIELD INCREMENTS INTO INCREMENTS AND PRODUCE THEM FOR THE SIZE-CLASSES IF NEEDED, ALSO SAVE INCREMENTS OF LOGS IN SEPARATE VARIABLES (THIS IS IN CASE ONLY TOTAL CHLOROPHYLL WAS ASSIMILATED)
60
61  IF (lwp) WRITE(numout,*) 'Code has gotten to point started'
62  CALL flush(numout)
63
64
65      IF ( ld_logchlpftinc ) THEN
66         !
67         ! Assimilating separate PFTs, so separately transform each from LogChl to Chl
68         !
69         IF ( npfts /= 4 ) THEN
70            CALL ctl_stop( 'If assimilating PFTs into ERSEM, nn_asmpfts must be 4' )
71         ENDIF
72         DO jj = 1, jpj
73            DO ji = 1, jpi
74               IF ( ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) > 0.0 ) .AND. &
75                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) > 0.0 ) .AND. &
76                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) > 0.0 ) .AND. &
77                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) > 0.0 ) ) THEN
78                  IF ( logchl_bkginc(ji,jj,1) /= 0.0 ) THEN
79                     logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) = 10**( LOG10( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) ) + &
80                        &                                                   logchl_bkginc(ji,jj,1) )                      - &
81                        &                                             trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1)
82                  ENDIF
83                  IF ( logchl_bkginc(ji,jj,2) /= 0.0 ) THEN
84                     logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) = 10**( LOG10( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) ) + &
85                        &                                                   logchl_bkginc(ji,jj,2) )                      - &
86                        &                                             trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2)
87                  ENDIF
88                  IF ( logchl_bkginc(ji,jj,3) /= 0.0 ) THEN
89                     logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) = 10**( LOG10( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) ) + &
90                        &                                                   logchl_bkginc(ji,jj,3) )                      - &
91                        &                                             trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3)
92                  ENDIF
93                  IF ( logchl_bkginc(ji,jj,4) /= 0.0 ) THEN
94                     logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) = 10**( LOG10( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) ) + &
95                        &                                                   logchl_bkginc(ji,jj,4) )                      - &
96                        &                                             trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4)
97                  ENDIF
98                  IF (k_maxchlinc > 0.0) THEN
99                     chl_inc = logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) + &
100                               logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) + &
101                               logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) + &
102                               logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl4)
103                     IF ( ABS(chl_inc) > k_maxchlinc ) THEN
104                        chl_tot = ABS(chl_inc) / k_maxchlinc
105                        logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) = logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) / chl_tot
106                        logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) = logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) / chl_tot
107                        logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) = logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) / chl_tot
108                        logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) = logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) / chl_tot
109                     ENDIF
110                  ENDIF
111               ENDIF
112            END DO
113         END DO
114      ELSE
115         !
116         ! Assimilating total Chl, so transform total from LogChl to Chl
117         ! and split between PFTs according to the existing background ratios
118         !
119         IF ( npfts /= 1 ) THEN
120            CALL ctl_stop( 'If assimilating total chlorophyll, nn_asmpfts must be 1' )
121         ENDIF
122         DO jj = 1, jpj
123            DO ji = 1, jpi
124               IF ( ( logchl_bkginc(ji,jj,1)               /= 0.0 ) .AND. &
125                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) >  0.0 ) .AND. &
126                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) >  0.0 ) .AND. &
127                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) >  0.0 ) .AND. &
128                  & ( trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) >  0.0 ) ) THEN
129                  chl_tot = trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) + trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) + &
130                     &      trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) + trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4)
131                  chl_inc = 10**( LOG10( chl_tot ) + logchl_bkginc(ji,jj,1) ) - chl_tot
132                  IF (k_maxchlinc > 0.0) THEN
133                     chl_inc = MAX( -1.0 * k_maxchlinc, MIN( chl_inc, k_maxchlinc ) )
134                  ENDIF
135                  logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) = chl_inc * trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) / chl_tot
136                  logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) = chl_inc * trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) / chl_tot
137                  logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) = chl_inc * trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) / chl_tot
138                  logchl_balinc(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) = chl_inc * trn(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) / chl_tot
139               ENDIF
140            END DO
141         END DO
142      ENDIF
143
144     
145     
146   
147
148
149
150  IF (lwp) WRITE(numout,*) 'Code has gotten to point log_inc transformed to inc'
151  CALL flush(numout)
152   ! STEP 2: READ CROSS_COVARIANCES
153
154  CALL read_crosscovariances(cc_chl_field, npfts)
155
156 
157 IF (lwp) WRITE(numout,*) 'Code has gotten to point cc read'
158 CALL flush(numout)
159  !! STEP 3: CALCULATE INCREMENTS OF OTHER BIOGEOCHEMICAL VARIABLES
160
161
162  threshold(jp_fabm_m1+jp_fabm_n3n) = 100
163  threshold(jp_fabm_m1+jp_fabm_n4n) = 10
164  threshold(jp_fabm_m1+jp_fabm_n1p) = 15
165  threshold(jp_fabm_m1+jp_fabm_n5s) = 80
166  threshold(jp_fabm_m1+jp_fabm_z4c) = 20
167  threshold(jp_fabm_m1+jp_fabm_z5c) = 20
168  threshold(jp_fabm_m1+jp_fabm_p1c) = 25
169  threshold(jp_fabm_m1+jp_fabm_p2c) = 25
170  threshold(jp_fabm_m1+jp_fabm_p3c) = 25
171  threshold(jp_fabm_m1+jp_fabm_p4c) = 25
172  threshold(jp_fabm_m1+jp_fabm_p1n) = 10
173  threshold(jp_fabm_m1+jp_fabm_p2n) = 10
174  threshold(jp_fabm_m1+jp_fabm_p3n) = 10
175  threshold(jp_fabm_m1+jp_fabm_p4n) = 10
176  threshold(jp_fabm_m1+jp_fabm_p1p) = 10
177  threshold(jp_fabm_m1+jp_fabm_p2p) = 10
178  threshold(jp_fabm_m1+jp_fabm_p3p) = 10
179  threshold(jp_fabm_m1+jp_fabm_p4p) = 10
180
181 
182
183 
184  DO jj = 1, jpj
185         DO ji = 1, jpi
186                 DO jp=1,jptra
187
188                 
189            !CALL t;ge'e;fdgfgd
190
191                         IF ( (jp /= jp_fabm_m1+jp_fabm_chl1) .AND. &
192                            & (jp /= jp_fabm_m1+jp_fabm_chl2) .AND. &
193                            & (jp /= jp_fabm_m1+jp_fabm_chl3) .AND. &
194                            & (jp /= jp_fabm_m1+jp_fabm_chl4)) THEN
195           
196                         CALL asm_inc_biogeo_ersem(npfts, jp, trn(ji,jj,1,jp), logchl_bkginc(ji,jj,:),  cc_chl_field(ji,jj,:,:), threshold(jp), &
197                          & logchl_balinc(ji,jj,1,jp))
198
199
200                         END IF
201
202
203         END DO
204
205      END DO
206
207    END DO
208
209
210
211 IF (lwp) WRITE(numout,*) 'Code has gotten to point surface logs calculated'
212
213
214
215
216 
217
218
219
220
221 CALL flush(numout)
222       !! STEP 4: PROPAGATE EACH INCREMENT THROUGH THE MIXED LAYER
223
224
225      SELECT CASE( mld_choice_bgc )
226      CASE ( 1 )                   ! Turbocline/mixing depth [W points]
227         zmld(:,:) = hmld(:,:)
228      CASE ( 2 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [W points]
229         zmld(:,:) = hmlp(:,:)
230      CASE ( 3 )                   ! Kara MLD [Interpolated]
231#if defined key_karaml
232         IF ( ln_kara ) THEN
233            zmld(:,:) = hmld_kara(:,:)
234         ELSE
235            CALL ctl_stop( ' Kara mixed layer requested for LogChl assimilation,', &
236               &           ' but ln_kara=.false.' )
237         ENDIF
238#else
239         CALL ctl_stop( ' Kara mixed layer requested for LogChl assimilation,', &
240            &           ' but is not defined' )
241#endif
242      CASE ( 4 )                   ! Temperature criterion (0.2 K change from surface) [T points]
243         zmld(:,:) = hmld_tref(:,:)
244      CASE ( 5 )                   ! Density criterion (0.01 kg/m^3 change from 10m) [T points]
245         zmld(:,:) = hmlpt(:,:)
246      END SELECT
247      !
248      ! Now set MLD to bottom of a level and propagate incs equally through mixed layer
249      !
250
251       DO jj = 1, jpj
252         DO ji = 1, jpi
253            !
254            jkmax = jpk-1
255            DO jk = jpk-1, 1, -1
256               IF ( ( zmld(ji,jj) >  gdepw_n(ji,jj,jk)   ) .AND. &
257                  & ( zmld(ji,jj) <= gdepw_n(ji,jj,jk+1) ) ) THEN
258                  zmld(ji,jj) = gdepw_n(ji,jj,jk+1)
259                  jkmax = jk
260               ENDIF
261            END DO
262            !
263            DO jk = 2, jkmax
264
265               DO jp = 1, jptra
266
267               !IF (lwp) WRITE(numout,*) logchl_balinc_ersem_chl1(ji,jj,1)
268               !CALL flush(numout)
269               logchl_balinc(ji,jj,jk,jp) = logchl_balinc(ji,jj,1,jp)
270             
271               END DO
272            END DO
273
274 
275         END DO
276      END DO
277
278
279   IF (lwp) WRITE(numout,*) 'Code has gotten to point log transformed into inc and propagated'
280   CALL flush(numout)
281
282END SUBROUTINE asm_logchl_bal_ersem
283
284
285
286!CONTAINS
287SUBROUTINE read_crosscovariances(cc_chl_field, npfts)
288
289  !!     status = nf90_open(path = filepath, cmode = nf90_nowrite, ncid = ccfile)
290
291 INTEGER, INTENT(in  ) :: npfts
292 REAL(wp), INTENT( out), DIMENSION(jpi,jpj,npfts,jptra) :: cc_chl_field  !! CROSS-COV BETWEEN SPECIFIC SIZE CLASS AND GIVEN BIOGEOCHEM VAR
293 INTEGER :: ncid
294 INTEGER :: pft
295 CHARACTER(LEN=256)  :: path_chl(npfts)
296 CHARACTER(LEN=256)  :: id_chl(npfts)
297
298 IF (npfts /= 4) THEN
299
300 IF (lwp) WRITE(numout,*) 'Wrong number of Pfts for reading crosscovs files'
301  CALL flush(numout)
302
303 ELSE
304
305 path_chl(1) = "crosscovs_chl1.nc"
306 path_chl(2) = "crosscovs_chl2.nc"
307 path_chl(3) = "crosscovs_chl3.nc"
308 path_chl(4) = "crosscovs_chl4.nc"
309
310 id_chl(1) = "cc_Chl1"
311 id_chl(2) = "cc_Chl2"
312 id_chl(3) = "cc_Chl3"
313 id_chl(4) = "cc_Chl4"
314
315
316 DO pft = 1, npfts
317
318 CALL iom_open( path_chl(pft), ncid )
319 CALL iom_get( ncid, jpdom_autoglo, 'cc_NO3', cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_n3n))
320 CALL iom_get( ncid, jpdom_autoglo, 'cc_NH4', cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_n4n))
321 CALL iom_get( ncid, jpdom_autoglo, "cc_PO4", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_n1p))
322 CALL iom_get( ncid, jpdom_autoglo, "cc_SiO", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_n5s))
323 CALL iom_get( ncid, jpdom_autoglo, "cc_Z4c", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_z4c))
324 CALL iom_get( ncid, jpdom_autoglo, "cc_Z5c", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_z5c))
325 CALL iom_get( ncid, jpdom_autoglo, "cc_P1c", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p1c))
326 CALL iom_get( ncid, jpdom_autoglo, "cc_P2c", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p2c))
327 CALL iom_get( ncid, jpdom_autoglo, "cc_P3c", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p3c))
328 CALL iom_get( ncid, jpdom_autoglo, "cc_P4c", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p4c))
329 CALL iom_get( ncid, jpdom_autoglo, "cc_P1n", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p1n))
330 CALL iom_get( ncid, jpdom_autoglo, "cc_P2n", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p2n))
331 CALL iom_get( ncid, jpdom_autoglo, "cc_P3n", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p3n))
332 CALL iom_get( ncid, jpdom_autoglo, "cc_P4n", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p4n))
333 CALL iom_get( ncid, jpdom_autoglo, "cc_P1p", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p1p))
334 CALL iom_get( ncid, jpdom_autoglo, "cc_P2p", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p2p))
335 CALL iom_get( ncid, jpdom_autoglo, "cc_P3p", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p3p))
336 CALL iom_get( ncid, jpdom_autoglo, "cc_P4p", cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_p4p))
337 CALL iom_get( ncid, jpdom_autoglo, id_chl(pft), cc_chl_field(:,:,pft,jp_fabm_m1+jp_fabm_chl1))
338 CALL iom_close( ncid )
339
340 END DO
341
342 END IF
343
344
345  END SUBROUTINE read_crosscovariances
346
347
348!CONTAINS
349 SUBROUTINE asm_inc_biogeo_ersem(npfts, jp, bgf_value, chl, cc, threshold, output)
350 INTEGER, INTENT(in ) :: npfts
351 INTEGER, INTENT(in ) :: jp
352 REAL, INTENT(in  )  :: bgf_value
353 REAL, INTENT(in   ), DIMENSION(npfts) :: chl
354 REAL, INTENT(in   ), DIMENSION(npfts,jptra) :: cc  !! VARIANCE OF CHL SIZE-CLASS
355 REAL, INTENT(in   ) :: threshold
356 REAL, INTENT(  out) :: output
357
358
359 REAL :: total, weight_1, weight_2, weight_3, weight_4
360
361
362    output = 0.0
363
364    IF ( cc(1,jp)*cc(2,jp)*cc(3,jp)*cc(4,jp)*cc(1, jp_fabm_m1+jp_fabm_chl1)*cc(2, jp_fabm_m1+jp_fabm_chl2)*cc(3, jp_fabm_m1+jp_fabm_chl3)*cc(4, jp_fabm_m1+jp_fabm_chl4) /= 0.0  ) THEN
365
366       total = ABS(cc(1,jp)/cc(1, jp_fabm_m1+jp_fabm_chl1)) + ABS(cc(2,jp)/cc(2, jp_fabm_m1+jp_fabm_chl2)) + ABS(cc(3,jp)/cc(3, jp_fabm_m1+jp_fabm_chl3)) + ABS(cc(4,jp)/cc(4, jp_fabm_m1+jp_fabm_chl4))
367       weight_1 = ABS(cc(1,jp)/cc(1, jp_fabm_m1+jp_fabm_chl1))/total
368       weight_2 = ABS(cc(2,jp)/cc(2, jp_fabm_m1+jp_fabm_chl2))/total
369       weight_3 = ABS(cc(3,jp)/cc(3, jp_fabm_m1+jp_fabm_chl3))/total
370       weight_4 = ABS(cc(4,jp)/cc(4, jp_fabm_m1+jp_fabm_chl4))/total
371
372!       PRINT *, 'weights', weight_1, weight_2, weight_3, weight_4
373     
374       output = weight_1*(cc(1,jp)/cc(1, jp_fabm_m1+jp_fabm_chl1))*chl(1) + weight_2*(cc(2,jp)/cc(2, jp_fabm_m1+jp_fabm_chl2))*chl(2) + weight_3*(cc(3,jp)/cc(3, jp_fabm_m1+jp_fabm_chl3))*chl(3) + weight_4*(cc(4,jp)/cc(4, jp_fabm_m1+jp_fabm_chl4))*chl(4)
375
376
377       ! FIRST BREAK = DON'T INCREASE VALUE MORE THAN TENFOLD
378       
379       output = MIN(output, LOG10(11.0))
380     
381
382       ! SECOND BREAK = IF VALUE LARGER THAN THRESHOLD DO NOT INCREASE IT
383
384       IF (bgf_value > threshold) THEN
385           output = MIN(output, 0.0)
386       ENDIF
387
388        ! THIRD BREAK = IF VALUE MORE THAN FIFTH OF THRESHOLD THAN DON'T INCREASE IT BY MORE THAN AN INCREMENT EQUAL TO TENTH OF THRESHOLD (INCREASE IT SLOWLY)
389
390       IF (bgf_value > 0.2*threshold) THEN
391           output = MIN(output, LOG10(0.1*threshold/bgf_value + 1))
392       ENDIF
393
394       output = (10**( output ) - 1)*bgf_value
395 
396       PRINT *, 'output', output
397         
398     
399    ENDIF
400
401     
402     !! calculate increments for one or all variables
403
404   END SUBROUTINE asm_inc_biogeo_ersem
405
406
407
408 
409
410
411#else
412   !!----------------------------------------------------------------------
413   !!   Default option : Empty routine
414   !!----------------------------------------------------------------------
415CONTAINS
416   CONTAINS
417   SUBROUTINE asm_logchl_bal_ersem( ld_logchlpftinc, npfts, mld_choice_bgc, &
418      &                             k_maxchlinc, logchl_bkginc, logchl_balinc )
419      LOGICAL :: ld_logchlpftinc
420      INTEGER :: npfts
421      INTEGER :: mld_choice_bgc
422      REAL    :: k_maxchlinc
423      REAL    :: logchl_bkginc(:,:,:)
424      REAL    :: logchl_balinc(:,:,:,:)
425      WRITE(*,*) 'asm_logchl_bal_ersem: You should not have seen this print! error?'
426   END SUBROUTINE asm_logchl_bal_ersem
427
428
429
430
431
432 
433
434#endif
435   !!======================================================================
436END MODULE asmlogchlbal_ersem
Note: See TracBrowser for help on using the repository browser.