source: LMDZ6/branches/blowing_snow/libf/phylmdiso/limit_read_mod.F90 @ 5213

Last change on this file since 5213 was 4143, checked in by dcugnet, 3 years ago
  • Some variables are renamed or replaced by direct equivalents:
    • iso_indnum -> tracers(:)%iso_iName
    • niso_possibles -> niso
    • iqiso -> iqIsoPha ; index_trac -> itZonIso
    • ok_iso_verif -> isoCheck
    • ntraceurs_zone -> nzone ; ntraciso -> ntiso
    • qperemin -> min_qparent ; masseqmin -> min_qmass ; ratiomin -> min_ratio
  • Some renamed variables are only aliased with the older name (using USE <module>, ONLY: <oldName> => <newName>) in routines where they are repeated many times.
  • Few hard-coded indexes are now computed (examples: ilic, iso, ivap, irneb, iq_vap, iq_liq, iso_H2O, iso_HDO, iso_HTO, iso_O17, iso_O18).
  • The IF(isoCheck) test is now embedded in the check_isotopes_seq and check_isotopes_loc routines (lighter calling).
File size: 22.2 KB
Line 
1!
2! $Id: limit_read_mod.F90 3435 2019-01-22 15:21:59Z fairhead $
3!
4MODULE limit_read_mod
5!
6! This module reads the fichier "limit.nc" containing fields for surface forcing.
7!
8! Module subroutines :
9!  limit_read_frac    : call limit_read_tot and return the fractions
10!  limit_read_rug_alb : return rugosity and albedo, if coupled ocean call limit_read_tot first
11!  limit_read_sst     : return sea ice temperature   
12!  limit_read_tot     : read limit.nc and store the fields in local modules variables
13!
14  IMPLICIT NONE
15
16  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE, PRIVATE :: pctsrf
17!$OMP THREADPRIVATE(pctsrf)
18  REAL, ALLOCATABLE, DIMENSION(:),   SAVE, PRIVATE :: rugos
19!$OMP THREADPRIVATE(rugos)
20  REAL, ALLOCATABLE, DIMENSION(:),   SAVE, PRIVATE :: albedo
21!$OMP THREADPRIVATE(albedo) 
22  REAL, ALLOCATABLE, DIMENSION(:),   SAVE, PRIVATE :: sst
23!$OMP THREADPRIVATE(sst)   
24#ifdef ISO
25  REAL, ALLOCATABLE, DIMENSION(:),   SAVE, PRIVATE :: tuoce
26!$OMP THREADPRIVATE(tuoce)
27#endif
28  LOGICAL,SAVE :: read_continents=.FALSE.
29!$OMP THREADPRIVATE(read_continents)
30
31CONTAINS
32!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
33!!
34!! Public subroutines :
35!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
36
37
38  SUBROUTINE init_limit_read(first_day)
39  USE mod_grid_phy_lmdz
40  USE surface_data
41  USE mod_phys_lmdz_para
42#ifdef CPP_XIOS
43  USE XIOS
44#endif
45  IMPLICIT NONE
46    INTEGER, INTENT(IN) :: first_day
47   
48   
49    IF ( type_ocean /= 'couple') THEN
50      IF (grid_type==unstructured) THEN
51#ifdef CPP_XIOS
52        IF (is_omp_master) CALL xios_set_file_attr("limit_read",enabled=.TRUE.,record_offset=first_day)
53#endif
54      ENDIF 
55    ENDIF
56
57  END SUBROUTINE init_limit_read
58 
59  SUBROUTINE limit_read_frac(itime, dtime, jour, pctsrf_new, is_modified)
60!
61! This subroutine is called from "change_srf_frac" for case of
62! ocean=force or from ocean_slab_frac for ocean=slab.
63! The fraction for all sub-surfaces at actual time step is returned.
64
65    USE dimphy
66    USE indice_sol_mod
67
68! Input arguments
69!****************************************************************************************
70    INTEGER, INTENT(IN) :: itime   ! time step
71    INTEGER, INTENT(IN) :: jour    ! current day
72    REAL   , INTENT(IN) :: dtime   ! length of time step
73 
74! Output arguments
75!****************************************************************************************
76    REAL, DIMENSION(klon,nbsrf), INTENT(OUT) :: pctsrf_new  ! sub surface fractions
77    LOGICAL, INTENT(OUT)                     :: is_modified ! true if pctsrf is modified at this time step
78
79! End declaration
80!****************************************************************************************
81
82! 1) Read file limit.nc
83    CALL limit_read_tot(itime, dtime, jour, is_modified)
84
85! 2) Return the fraction read in limit_read_tot
86    pctsrf_new(:,:) = pctsrf(:,:)
87   
88  END SUBROUTINE limit_read_frac
89
90!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
91
92  SUBROUTINE limit_read_rug_alb(itime, dtime, jour, &
93       knon, knindex, &
94       rugos_out, alb_out)
95!
96! This subroutine is called from surf_land_bucket.
97! The flag "ok_veget" must can not be true. If coupled run, "ocean=couple"
98! then this routine will call limit_read_tot.
99!
100    USE dimphy
101    USE surface_data
102#ifdef ISO
103    USE isotopes_mod, ONLY: P_veg
104#endif
105
106! Input arguments
107!****************************************************************************************
108    INTEGER, INTENT(IN) :: itime                     ! numero du pas de temps courant
109    INTEGER, INTENT(IN) :: jour                      ! jour a lire dans l'annee
110    REAL   , INTENT(IN) :: dtime                     ! pas de temps de la physique (en s)
111    INTEGER, INTENT(IN) :: knon                      ! nomber of points on compressed grid
112    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex  ! grid point number for compressed grid
113! Output arguments
114!****************************************************************************************
115    REAL, DIMENSION(klon), INTENT(OUT) :: rugos_out
116    REAL, DIMENSION(klon), INTENT(OUT) :: alb_out
117   
118! Local variables
119!****************************************************************************************
120    INTEGER :: i
121    LOGICAL :: is_modified
122
123!****************************************************************************************
124
125IF (type_ocean == 'couple'.OR. &
126         (type_ocean == 'slab' .AND. version_ocean == 'sicINT')) THEN
127       ! limit.nc has not yet been read. Do it now!
128       CALL limit_read_tot(itime, dtime, jour, is_modified)
129    END IF
130
131    DO i=1,knon
132       rugos_out(i) = rugos(knindex(i))
133       alb_out(i)  = albedo(knindex(i))
134    END DO
135
136  END SUBROUTINE limit_read_rug_alb
137
138!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
139
140  SUBROUTINE limit_read_sst(knon, knindex, sst_out &
141#ifdef ISO
142     &  ,Roce,rlat   &
143#endif           
144    )
145!
146! This subroutine returns the sea surface temperature already read from limit.nc.
147!
148    USE dimphy, ONLY : klon
149#ifdef ISO
150    USE infotrac_phy, ONLY: niso
151    USE isotopes_mod, ONLY: tcorr,toce,modif_sst, &
152   &    deltaTtest,sstlatcrit,deltaTtestpoles,dsstlatcrit, &
153   &    iso_HTO,ok_prod_nucl_tritium
154#ifdef ISOVERIF
155    USE isotopes_verif_mod, ONLY: iso_verif_egalite_vect2D,iso_verif_positif, &
156        iso_verif_positif_nostop
157#endif
158#endif
159
160    INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
161    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex  ! grid point number for compressed grid
162    REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out
163#ifdef ISO
164    REAL, DIMENSION(klon) :: tuoce_out ! sortie tritium surface ocean
165#endif
166
167    INTEGER :: i
168#ifdef ISO
169  real, intent(out), dimension(niso,klon) :: Roce
170  integer :: ixt 
171  real, intent(in),dimension(klon) :: rlat
172  real lat_locale
173#endif
174!#ifdef ISOVERIF
175!   integer iso_verif_positif_nostop
176!#endif
177
178    DO i = 1, knon
179       sst_out(i) = sst(knindex(i))
180    END DO
181
182
183#ifdef ISO
184     if (iso_HTO.gt.0) then
185     if (ok_prod_nucl_tritium) then ! si on active la production nucleaire de tritium
186        DO i = 1, knon
187          tuoce_out(i)=tuoce(knindex(i))
188        END DO
189     endif
190     endif
191#endif
192
193#ifdef ISO
194  if (modif_sst.ge.1) then
195  do i = 1, knon
196    lat_locale=rlat(knindex(i)) 
197    ! test: modification uniforme de la sst
198    if (modif_sst.eq.1) then   
199       sst_out(i)= sst_out(i)+deltaTtest 
200    elseif (modif_sst.eq.2) then   !if (modif_sst.eq.1) then       
201        ! pattern parabolique en dehors des tropiques (sstlatcrit)
202        if (abs(lat_locale).gt.sstlatcrit) then
203          sst_out(i)= sst_out(i)+deltaTtestpoles &
204   &             *(lat_locale**2-sstlatcrit**2) &
205   &             /(90.0**2-sstlatcrit**2)                     
206        endif !if (abs(lat_locale).gt.abs(sstlatcrit)) then
207
208    else if (modif_sst.eq.3) then
209
210        if (abs(lat_locale).gt.abs(sstlatcrit)) then
211            if (abs(lat_locale).gt.sstlatcrit+dsstlatcrit) then
212                sst_out(i)= sst_out(i)+deltaTtestpoles
213            else
214                sst_out(i)= sst_out(i)+deltaTtestpoles &
215    &               *(abs(lat_locale)-sstlatcrit)/dsstlatcrit
216            endif
217        endif
218    endif !if (modif_sst.eq.1) then   
219    enddo !do i = 1, knon
220    endif !if (modif_sst.ge.1) then
221#endif
222#ifdef ISOVERIF
223     do i=1,knon,20
224       call iso_verif_positif(sst_out(i)-100.0,'limit_read 4323')
225     enddo
226#endif
227
228
229#ifdef ISO
230        !* lecture de Roce
231        ! 1) première possibilité: valeur fixe à SMOW
232        DO i = 1, knon
233          do ixt=1,niso
234            Roce(ixt,i)=tcorr(ixt)*toce(ixt)
235          enddo !do ixt=1,niso
236        enddo !DO i = 1, knon
237        ! 2) deuxième possibilité: lecture de la carte
238        ! A FAIRE
239
240        ! lecture pour le tritium
241        if ((iso_HTO.gt.0).and.(ok_prod_nucl_tritium)) then
242             ! lecture de la carte tritium ocean surface
243             Roce(iso_HTO,i)=tcorr(iso_HTO)*tuoce_out(i)*1.E-18*2.
244        endif       
245#endif 
246
247#ifdef ISOVERIF
248        do i=1,knon
249          if (iso_verif_positif_nostop(370.0-sst_out(i), &
250              'limit_read 368').eq.1) then
251             write(*,*) 'i,knindex,sst_out=',i,knindex,sst_out(i)
252             stop
253          endif
254        enddo
255#endif     
256
257
258  END SUBROUTINE limit_read_sst
259
260!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
261!!
262!! Private subroutine :
263!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
264
265  SUBROUTINE limit_read_tot(itime, dtime, jour, is_modified)
266!
267! Read everything needed from limit.nc
268!
269! 0) Initialize
270! 1) Open the file limit.nc, if it is time
271! 2) Read fraction, if not type_ocean=couple
272! 3) Read sea surface temperature, if not type_ocean=couple
273! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
274! 5) Close file and distribuate variables to all processus
275
276    USE dimphy
277    USE mod_grid_phy_lmdz
278    USE mod_phys_lmdz_para
279    USE surface_data, ONLY : type_ocean, ok_veget
280    USE netcdf
281    USE indice_sol_mod
282#ifdef ISO
283    USE isotopes_mod, ONLY : iso_HTO,ok_prod_nucl_tritium
284#ifdef ISOVERIF
285    USE isotopes_verif_mod, ONLY : iso_verif_positif_nostop
286#endif
287#endif
288    USE phys_cal_mod, ONLY : calend, year_len
289    USE print_control_mod, ONLY: lunout, prt_level
290#ifdef CPP_XIOS
291    USE XIOS, ONLY: xios_recv_field
292#endif
293   
294    IMPLICIT NONE
295   
296! In- and ouput arguments
297!****************************************************************************************
298    INTEGER, INTENT(IN) :: itime   ! numero du pas de temps courant
299    INTEGER, INTENT(IN) :: jour    ! jour a lire dans l'annee
300    REAL   , INTENT(IN) :: dtime   ! pas de temps de la physique (en s)
301
302    LOGICAL, INTENT(OUT) :: is_modified  ! true if pctsrf is modified at this time step
303
304! Locals variables with attribute SAVE
305!****************************************************************************************
306! frequence de lecture des conditions limites (en pas de physique)
307    INTEGER,SAVE                              :: lmt_pas
308!$OMP THREADPRIVATE(lmt_pas)
309    LOGICAL, SAVE                             :: first_call=.TRUE.
310!$OMP THREADPRIVATE(first_call) 
311    INTEGER, SAVE                             :: jour_lu = -1
312!$OMP THREADPRIVATE(jour_lu) 
313! Locals variables
314!****************************************************************************************
315    INTEGER                                   :: nid, nvarid, ndimid, nn
316    INTEGER                                   :: ii, ierr
317    INTEGER, DIMENSION(2)                     :: start, epais
318    REAL, DIMENSION(klon_glo,nbsrf)           :: pct_glo  ! fraction at global grid
319    REAL, DIMENSION(klon_glo)                 :: sst_glo  ! sea-surface temperature at global grid
320    REAL, DIMENSION(klon_glo)                 :: rug_glo  ! rugosity at global grid
321    REAL, DIMENSION(klon_glo)                 :: alb_glo  ! albedo at global grid
322
323    REAL, DIMENSION(klon_mpi,nbsrf)           :: pct_mpi  ! fraction at global grid
324    REAL, DIMENSION(klon_mpi)                 :: sst_mpi  ! sea-surface temperature at global grid
325    REAL, DIMENSION(klon_mpi)                 :: rug_mpi  ! rugosity at global grid
326    REAL, DIMENSION(klon_mpi)                 :: alb_mpi  ! albedo at global grid
327
328    CHARACTER(len=20)                         :: modname='limit_read_mod'     
329    CHARACTER(LEN=99)                         :: abort_message, calendar, str
330#ifdef ISO
331    REAL, DIMENSION(klon_glo)                 :: tuoce_glo  ! sea-surface tritium et global grid
332#endif
333
334! End declaration
335!****************************************************************************************
336
337!****************************************************************************************
338! 0) Initialization
339!
340!****************************************************************************************
341    IF (first_call) THEN
342       first_call=.FALSE.
343       ! calculate number of time steps for one day
344       lmt_pas = NINT(86400./dtime * 1.0)
345       
346       ! Allocate module save variables
347       IF ( type_ocean /= 'couple' ) THEN
348          ALLOCATE(pctsrf(klon,nbsrf), sst(klon), stat=ierr)
349          IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating pctsrf and sst',1)
350       END IF
351
352       IF ( .NOT. ok_veget ) THEN
353          ALLOCATE(rugos(klon), albedo(klon), stat=ierr)
354          IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating rugos and albedo',1)
355       END IF
356
357!$OMP MASTER  ! Only master thread
358       IF (is_mpi_root) THEN ! Only master processus
359          ierr = NF90_OPEN ('limit.nc', NF90_NOWRITE, nid)
360          IF (ierr /= NF90_NOERR) CALL abort_physic(modname,&
361               'Pb d''ouverture du fichier de conditions aux limites',1)
362
363          !--- WARNING IF CALENDAR IS KNOWN AND DOES NOT MATCH THE ONE OF LMDZ
364          ierr=NF90_INQ_VARID(nid, 'TEMPS', nvarid)
365          ierr=NF90_GET_ATT(nid, nvarid, 'calendar', calendar)
366          IF(ierr==NF90_NOERR.AND.calendar/=calend.AND.prt_level>=1) THEN
367             WRITE(lunout,*)'BEWARE: gcm and limit.nc calendars differ: '
368             WRITE(lunout,*)'  '//TRIM(calend)//' for gcm'
369             WRITE(lunout,*)'  '//TRIM(calendar)//' for limit.nc file'
370          END IF
371
372          !--- ERROR IF FILE RECORDS NUMBER IS NOT EQUAL TO EXPECTED NUMBER OF DAYS         
373          IF (grid_type==unstructured) THEN
374            ierr=NF90_INQ_DIMID(nid,"time_year",ndimid)
375          ELSE
376            ierr=NF90_INQUIRE(nid, UnlimitedDimID=ndimid)
377          ENDIF
378          ierr=NF90_INQUIRE_DIMENSION(nid, ndimid, len=nn)
379          WRITE(abort_message,'(a,2(i0,a))')'limit.nc records number (',nn,') does no'//&
380            't match year length (',year_len,')'
381          IF(nn/=year_len) CALL abort_physic(modname,abort_message,1)
382
383          !--- ERROR IF FILES AND LMDZ HORIZONTAL RESOLUTIONS DO NOT MATCH
384          IF (grid_type==unstructured) THEN
385            ierr=NF90_INQ_DIMID(nid, 'cell', ndimid)
386          ELSE
387            ierr=NF90_INQ_DIMID(nid, 'points_physiques', ndimid)
388          ENDIF
389          ierr=NF90_INQUIRE_DIMENSION(nid, ndimid, len=nn)
390          WRITE(abort_message,'(a,2(i0,a))')'limit.nc horizontal number of cells (',nn, &
391            ') does not match LMDZ klon_glo (',klon_glo,')'
392          IF(nn/=klon_glo) CALL abort_physic(modname,abort_message,1)
393
394          ierr = NF90_CLOSE(nid)
395          IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Pb when closing file', 1)
396       END IF ! is_mpi_root
397!$OMP END MASTER
398!$OMP BARRIER
399    END IF
400
401!****************************************************************************************
402! 1) Open the file limit.nc if it is the right moment to read, once a day.
403!    The file is read only by the master thread of the master mpi process(is_mpi_root)
404!    Check by the way if the number of records is correct.
405!
406!****************************************************************************************
407
408    is_modified = .FALSE.
409!ym    IF (MOD(itime-1, lmt_pas) == 0 .OR. jour_lu /= jour ) THEN   ! time to read
410!  not REALLY PERIODIC
411    IF (MOD(itime-1, lmt_pas) == 0 .OR. (jour_lu /= jour .AND. grid_type /= unstructured)) THEN   ! time to read
412!    IF (MOD(itime-1, lmt_pas) == 0) THEN   ! time to read
413       jour_lu = jour
414       is_modified = .TRUE.
415
416      IF (grid_type==unstructured) THEN
417
418#ifdef CPP_XIOS
419        IF ( type_ocean /= 'couple') THEN
420
421           IF (is_omp_master) CALL xios_recv_field("foce_limin",pct_mpi(:,is_oce))
422           IF (is_omp_master) CALL xios_recv_field("fsic_limin",pct_mpi(:,is_sic))
423  !         IF (read_continents .OR. itime == 1) THEN
424           IF (is_omp_master) CALL xios_recv_field("fter_limin",pct_mpi(:,is_ter))
425           IF (is_omp_master) CALL xios_recv_field("flic_limin",pct_mpi(:,is_lic))
426  !         ENDIF
427         ENDIF! type_ocean /= couple
428         
429         IF ( type_ocean /= 'couple') THEN                   
430             IF (is_omp_master) CALL xios_recv_field("sst_limin",sst_mpi)
431         ENDIF
432       
433         IF (.NOT. ok_veget) THEN
434           IF (is_omp_master) CALL xios_recv_field("alb_limin",alb_mpi)
435           IF (is_omp_master) CALL xios_recv_field("rug_limin",rug_mpi)
436         ENDIF
437
438       IF ( type_ocean /= 'couple') THEN
439          CALL Scatter_omp(sst_mpi,sst)
440          CALL Scatter_omp(pct_mpi(:,is_oce),pctsrf(:,is_oce))
441          CALL Scatter_omp(pct_mpi(:,is_sic),pctsrf(:,is_sic))
442!          IF (read_continents .OR. itime == 1) THEN
443             CALL Scatter_omp(pct_mpi(:,is_ter),pctsrf(:,is_ter))
444             CALL Scatter_omp(pct_mpi(:,is_lic),pctsrf(:,is_lic))
445!          END IF
446       END IF
447
448       IF (.NOT. ok_veget) THEN
449          CALL Scatter_omp(alb_mpi, albedo)
450          CALL Scatter_omp(rug_mpi, rugos)
451       END IF
452#endif
453
454 
455     ELSE      ! grid_type==regular
456
457!$OMP MASTER  ! Only master thread
458       IF (is_mpi_root) THEN ! Only master processus!
459
460          ierr = NF90_OPEN ('limit.nc', NF90_NOWRITE, nid)
461          IF (ierr /= NF90_NOERR) CALL abort_physic(modname,&
462               'Pb d''ouverture du fichier de conditions aux limites',1)
463
464          ! La tranche de donnees a lire:
465          start(1) = 1
466          start(2) = jour
467          epais(1) = klon_glo
468          epais(2) = 1
469
470
471!****************************************************************************************
472! 2) Read fraction if not type_ocean=couple
473!
474!****************************************************************************************
475
476          IF ( type_ocean /= 'couple') THEN
477!
478! Ocean fraction
479             ierr = NF90_INQ_VARID(nid, 'FOCE', nvarid)
480             IF (ierr /= NF90_NOERR) CALL abort_physic(modname, 'Le champ <FOCE> est absent',1)
481             
482             ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_oce),start,epais)
483             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FOCE>' ,1)
484!
485! Sea-ice fraction
486             ierr = NF90_INQ_VARID(nid, 'FSIC', nvarid)
487             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FSIC> est absent',1)
488
489             ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_sic),start,epais)
490             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FSIC>' ,1)
491
492
493! Read land and continentals fraction only if asked for
494             IF (read_continents .OR. itime == 1) THEN
495!
496! Land fraction
497                ierr = NF90_INQ_VARID(nid, 'FTER', nvarid)
498                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FTER> est absent',1)
499               
500                ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_ter),start,epais)
501                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FTER>',1)
502!
503! Continentale ice fraction
504                ierr = NF90_INQ_VARID(nid, 'FLIC', nvarid)
505                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FLIC> est absent',1)
506
507                ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_lic),start,epais)
508                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FLIC>',1)
509             END IF
510
511          END IF ! type_ocean /= couple
512
513!****************************************************************************************
514! 3) Read sea-surface temperature, if not coupled ocean
515!
516!****************************************************************************************
517          IF ( type_ocean /= 'couple') THEN
518
519             ierr = NF90_INQ_VARID(nid, 'SST', nvarid)
520             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <SST> est absent',1)
521
522             ierr = NF90_GET_VAR(nid,nvarid,sst_glo,start,epais)
523             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <SST>',1)
524         
525#ifdef ISO
526             IF ((iso_HTO.gt.0).and.(ok_prod_nucl_tritium)) THEN
527               ierr = NF90_INQ_VARID(nid, 'TUOCE', nvarid)
528               IF (ierr /= NF90_NOERR) CALL abort_gcm(modname,'Le champ <TUOCE> est absent',1)
529
530               ierr = NF90_GET_VAR(nid,nvarid,tuoce_glo,start,epais)
531               IF (ierr /= NF90_NOERR) CALL abort_gcm(modname,'Lecture echouee pour <TUOCE>',1)
532             END IF
533#ifdef ISOVERIF
534             do ii=1,klon_glo
535               if (iso_verif_positif_nostop(370.0-sst_glo(ii),  &
536                'limit_read 384').eq.1) then
537                 write(*,*) 'ii,sst_glo=',ii,sst_glo(ii)
538                 write(*,*) 'jour,start,epais=',jour,start,epais
539                 stop
540               endif
541             enddo
542#endif
543#endif
544
545          END IF
546
547!****************************************************************************************
548! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
549!
550!****************************************************************************************
551
552          IF (.NOT. ok_veget) THEN
553!
554! Read albedo
555             ierr = NF90_INQ_VARID(nid, 'ALB', nvarid)
556             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <ALB> est absent',1)
557
558             ierr = NF90_GET_VAR(nid,nvarid,alb_glo,start,epais)
559             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <ALB>',1)
560!
561! Read rugosity
562             ierr = NF90_INQ_VARID(nid, 'RUG', nvarid)
563             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <RUG> est absent',1)
564
565             ierr = NF90_GET_VAR(nid,nvarid,rug_glo,start,epais)
566             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <RUG>',1)
567
568          END IF
569
570!****************************************************************************************
571! 5) Close file and distribuate variables to all processus
572!
573!****************************************************************************************
574          ierr = NF90_CLOSE(nid)
575          IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Pb when closing file', 1)
576       ENDIF ! is_mpi_root
577
578!$OMP END MASTER
579!$OMP BARRIER
580
581       IF ( type_ocean /= 'couple') THEN
582          CALL Scatter(sst_glo,sst)
583          CALL Scatter(pct_glo(:,is_oce),pctsrf(:,is_oce))
584          CALL Scatter(pct_glo(:,is_sic),pctsrf(:,is_sic))
585          IF (read_continents .OR. itime == 1) THEN
586             CALL Scatter(pct_glo(:,is_ter),pctsrf(:,is_ter))
587             CALL Scatter(pct_glo(:,is_lic),pctsrf(:,is_lic))
588          END IF
589#ifdef ISO
590          IF ((iso_HTO.gt.0).and.(ok_prod_nucl_tritium)) THEN
591             CALL Scatter(tuoce_glo,tuoce)
592          END IF
593#endif
594       END IF
595
596       IF (.NOT. ok_veget) THEN
597          CALL Scatter(alb_glo, albedo)
598          CALL Scatter(rug_glo, rugos)
599       END IF
600
601      ENDIF ! Grid type
602
603    ENDIF ! time to read
604
605  END SUBROUTINE limit_read_tot
606
607END MODULE limit_read_mod
Note: See TracBrowser for help on using the repository browser.