source: LMDZ.3.3/branches/rel-LF/libf/dyn3d/limit_netcdf.F @ 359

Last change on this file since 359 was 359, checked in by lmdzadmin, 22 years ago

Menage sur les etats initiaux
LF

  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 34.1 KB
Line 
1C
2C $Header$
3C
4      SUBROUTINE limit_netcdf(interbar, extrap, oldice, masque, pctsrf)
5c
6      IMPLICIT none
7c
8c-------------------------------------------------------------
9C Author : L. Fairhead
10C Date   : 27/01/94
11C Objet  : Construction des fichiers de conditions aux limites
12C          pour le nouveau
13C          modele a partir de fichiers de climatologie. Les deux
14C          grilles doivent etre regulieres
15c
16c Modifie par z.x.li (le23mars1994)
17c Modifie par L. Fairhead (fairhead@lmd.jussieu.fr) septembre 1999
18c                         pour lecture netcdf dans LMDZ.3.3
19c Modifie par P;Le Van  ,  juillet 2001
20c-------------------------------------------------------------
21c
22#include "dimensions.h"
23#include "paramet.h"
24#include "control.h"
25#include "logic.h"
26#include "netcdf.inc"
27#include "comvert.h"
28#include "comgeom2.h"
29#include "comconst.h"
30#include "dimphy.h"
31#include "indicesol.h"
32c
33c-----------------------------------------------------------------------
34      LOGICAL interbar, extrap, oldice
35
36      REAL phy_nat(klon,360), phy_nat0(klon)
37      REAL phy_alb(klon,360)
38      REAL phy_sst(klon,360)
39      REAL phy_bil(klon,360)
40      REAL phy_rug(klon,360)
41      REAL phy_ice(klon)
42c
43      real pctsrf_t(klon,nbsrf,360)
44      real pctsrf(klon,nbsrf)
45      REAL verif
46
47      REAL masque(iip1,jjp1)
48      REAL mask(iim,jjp1)
49CPB
50C newlmt indique l'utilisation de la sous-maille fractionnelle
51C tandis que l'ancien codage utilise l'indicateur du sol (0,1,2,3)
52      LOGICAL newlmt, fracterre
53      PARAMETER(newlmt=.TRUE.)
54      PARAMETER(fracterre = .TRUE.)
55
56C Declarations pour le champ de depart
57      INTEGER imdep, jmdep,lmdep
58      INTEGER  tbid
59      PARAMETER ( tbid = 60 )        ! >52 semaines
60      REAL  timecoord(tbid)
61c
62      REAL , ALLOCATABLE :: dlon_msk(:), dlat_msk(:)
63      REAL , ALLOCATABLE :: lonmsk_ini(:), latmsk_ini(:)
64      REAL , ALLOCATABLE :: dlon(:), dlat(:)
65      REAL , ALLOCATABLE :: dlon_ini(:), dlat_ini(:)
66      REAL , ALLOCATABLE :: champ_msk(:), champ(:)
67      REAL , ALLOCATABLE :: work(:,:)
68
69      CHARACTER*25 title
70
71C Declarations pour le champ interpole 2D
72      REAL champint(iim,jjp1)
73      real chmin,chmax
74
75C Declarations pour le champ interpole 3D
76      REAL champtime(iim,jjp1,tbid)
77      REAL timeyear(tbid)
78      REAL champan(iip1,jjp1,366)
79
80C Declarations pour l'inteprolation verticale
81      REAL ax(tbid), ay(tbid)
82      REAL by
83      REAL yder(tbid)
84
85
86      INTEGER ierr
87      INTEGER dimfirst(3)
88      INTEGER dimlast(3)
89c
90      INTEGER nid, ndim, ntim
91      INTEGER dims(2), debut(2), epais(2)
92      INTEGER id_tim
93      INTEGER id_NAT, id_SST, id_BILS, id_RUG, id_ALB
94CPB
95      INTEGER id_FOCE, id_FSIC, id_FTER, id_FLIC
96
97      INTEGER i, j, k, l, ji
98c declarations pour lecture glace de mer
99      INTEGER :: iml_lic, jml_lic, llm_tmp, ttm_tmp, iret
100      INTEGER :: itaul(1), fid
101      REAL :: lev(1), date, dt
102      REAL, ALLOCATABLE, DIMENSION(:,:) :: lon_lic, lat_lic
103      REAL, ALLOCATABLE, DIMENSION(:)  :: dlon_lic, dlat_lic
104      REAL, ALLOCATABLE, DIMENSION (:,:) :: fraclic
105      REAL :: flic_tmp(iip1, jjp1)
106
107c Diverses variables locales
108      REAL time
109! pour la lecture du fichier masque ocean
110      integer :: nid_o2a
111      logical :: couple = .false.
112      INTEGER :: iml_omask, jml_omask
113      REAL, ALLOCATABLE, DIMENSION(:,:) :: lon_omask, lat_omask
114      REAL, ALLOCATABLE, DIMENSION(:)  :: dlon_omask, dlat_omask
115      REAL, ALLOCATABLE, DIMENSION (:,:) :: ocemask, ocetmp
116      real, dimension(klon) :: ocemask_fi
117
118      INTEGER          longcles
119      PARAMETER      ( longcles = 20 )
120      REAL  clesphy0 ( longcles      )
121#include "serre.h"
122      INTEGER ncid,varid,ndimid(4),dimid
123      character*30 namedim
124      CHARACTER*80 :: varname
125
126c initialisations:
127      CALL conf_gcm( 99, .TRUE. , clesphy0 )
128
129
130      pi     = 4. * ATAN(1.)
131      rad    = 6 371 229.
132      omeg   = 4.* ASIN(1.)/(24.*3600.)
133      g      = 9.8
134      daysec = 86400.
135      kappa  = 0.2857143
136      cpp    = 1004.70885
137      dtvr    = daysec/FLOAT(day_step)
138      CALL inigeom
139c
140C Traitement du relief au sol
141c
142      write(*,*) 'Traitement du relief au sol pour fabriquer masque'
143      ierr = NF_OPEN('Relief.nc', NF_NOWRITE, ncid)
144
145      if (ierr.ne.0) then
146        print *, NF_STRERROR(ierr)
147        STOP
148      ENDIF
149
150      ierr = NF_INQ_VARID(ncid,'RELIEF',varid)
151      if (ierr.ne.0) then
152        print *, NF_STRERROR(ierr)
153        STOP
154      ENDIF
155      ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
156      if (ierr.ne.0) then
157        print *, NF_STRERROR(ierr)
158        STOP
159      ENDIF
160      ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
161      if (ierr.ne.0) then
162        print *, NF_STRERROR(ierr)
163        STOP
164      ENDIF
165      print*,'variable ', namedim, 'dimension ', imdep
166      ierr = NF_INQ_VARID(ncid,namedim,dimid)
167      if (ierr.ne.0) then
168        print *, NF_STRERROR(ierr)
169        STOP
170      ENDIF
171
172      ALLOCATE( lonmsk_ini(imdep) )
173      ALLOCATE(   dlon_msk(imdep) )
174
175#ifdef NC_DOUBLE
176      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,lonmsk_ini)
177#else
178      ierr = NF_GET_VAR_REAL(ncid,dimid,lonmsk_ini)
179#endif
180
181c
182      if (ierr.ne.0) then
183        print *, NF_STRERROR(ierr)
184        STOP
185      ENDIF
186      ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
187      if (ierr.ne.0) then
188        print *, NF_STRERROR(ierr)
189        STOP
190      ENDIF
191      print*,'variable ', namedim, 'dimension ', jmdep
192      ierr = NF_INQ_VARID(ncid,namedim,dimid)
193      if (ierr.ne.0) then
194        print *, NF_STRERROR(ierr)
195        STOP
196      ENDIF
197
198      ALLOCATE( latmsk_ini(jmdep) )
199      ALLOCATE(   dlat_msk(jmdep) )
200      ALLOCATE(  champ_msk(imdep*jmdep) )
201
202#ifdef NC_DOUBLE
203      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,latmsk_ini)
204#else
205      ierr = NF_GET_VAR_REAL(ncid,dimid,latmsk_ini)
206#endif
207c
208      if (ierr.ne.0) then
209        print *, NF_STRERROR(ierr)
210        STOP
211      ENDIF
212#ifdef NC_DOUBLE
213      ierr = NF_GET_VAR_DOUBLE(ncid,varid,champ_msk)
214#else
215      ierr = NF_GET_VAR_REAL(ncid,varid,champ_msk)
216#endif
217c
218      if (ierr.ne.0) then
219        print *, NF_STRERROR(ierr)
220        STOP
221      ENDIF
222c
223      title='RELIEF'
224
225      CALL conf_dat2d(title,imdep, jmdep, lonmsk_ini, latmsk_ini,
226     . dlon_msk, dlat_msk, champ_msk, interbar  )
227
228      DO i = 1, iim
229      DO j = 1, jjp1
230         mask(i,j) = masque(i,j)
231      ENDDO
232      ENDDO
233      WRITE(*,*) 'MASK:'
234      WRITE(*,'(96i1)')INT(mask)     
235      ierr = NF_CLOSE(ncid)
236c
237c
238C Traitement de la rugosite
239c
240      PRINT*, 'Traitement de la rugosite'
241      ierr = NF_OPEN('Rugos.nc', NF_NOWRITE, ncid)
242      if (ierr.ne.0) then
243        print *, NF_STRERROR(ierr)
244        STOP
245      ENDIF
246
247      ierr = NF_INQ_VARID(ncid,'RUGOS',varid)
248      if (ierr.ne.0) then
249        print *, NF_STRERROR(ierr)
250        STOP
251      ENDIF
252      ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
253      if (ierr.ne.0) then
254        print *, NF_STRERROR(ierr)
255        STOP
256      ENDIF
257      ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
258      if (ierr.ne.0) then
259        print *, NF_STRERROR(ierr)
260        STOP
261      ENDIF
262      print*,'variable ', namedim, 'dimension ', imdep
263      ierr = NF_INQ_VARID(ncid,namedim,dimid)
264      if (ierr.ne.0) then
265        print *, NF_STRERROR(ierr)
266        STOP
267      ENDIF
268
269      ALLOCATE( dlon_ini(imdep) )
270      ALLOCATE(     dlon(imdep) )
271
272#ifdef NC_DOUBLE
273      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
274#else
275      ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
276#endif
277      if (ierr.ne.0) then
278        print *, NF_STRERROR(ierr)
279        STOP
280      ENDIF
281      ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
282      if (ierr.ne.0) then
283        print *, NF_STRERROR(ierr)
284        STOP
285      ENDIF
286      print*,'variable ', namedim, 'dimension ', jmdep
287      ierr = NF_INQ_VARID(ncid,namedim,dimid)
288      if (ierr.ne.0) then
289        print *, NF_STRERROR(ierr)
290        STOP
291      ENDIF
292
293      ALLOCATE( dlat_ini(jmdep) )
294      ALLOCATE(     dlat(jmdep) )
295
296#ifdef NC_DOUBLE
297      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
298#else
299      ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
300#endif
301      if (ierr.ne.0) then
302        print *, NF_STRERROR(ierr)
303        STOP
304      ENDIF
305      ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
306      if (ierr.ne.0) then
307        print *, NF_STRERROR(ierr)
308        STOP
309      ENDIF
310      print*,'variable ', namedim, 'dimension ', lmdep
311      ierr = NF_INQ_VARID(ncid,namedim,dimid)
312      if (ierr.ne.0) then
313        print *, NF_STRERROR(ierr)
314        STOP
315      ENDIF
316#ifdef NC_DOUBLE
317      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
318#else
319      ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
320#endif
321      if (ierr.ne.0) then
322        print *, NF_STRERROR(ierr)
323        STOP
324      ENDIF
325c
326      ALLOCATE( champ(imdep*jmdep) )
327
328      DO  200 l = 1, lmdep
329         dimfirst(1) = 1
330         dimfirst(2) = 1
331         dimfirst(3) = l
332c
333         dimlast(1) = imdep
334         dimlast(2) = jmdep
335         dimlast(3) = 1
336c
337         PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
338         print*,dimfirst,dimlast
339#ifdef NC_DOUBLE
340         ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
341#else
342         ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
343#endif
344         if (ierr.ne.0) then
345           print *, NF_STRERROR(ierr)
346           STOP
347         ENDIF
348   
349        title = 'Rugosite Amip '
350c
351        CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
352     .                      dlon, dlat, champ, interbar          )
353
354       IF ( interbar )   THEN
355         DO j = 1, imdep * jmdep
356           champ(j) = LOG(champ(j))
357         ENDDO
358
359        IF( l.EQ.1 )  THEN
360         WRITE(6,*) '-------------------------------------------------',
361     ,'------------------------'
362         WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
363     , ' pour la rugosite $$$ '
364         WRITE(6,*) '-------------------------------------------------',
365     ,'------------------------'
366        ENDIF
367        CALL inter_barxy ( imdep,jmdep -1,dlon,dlat,champ ,
368     ,                  iim,jjm,rlonu,rlatv, jjp1,champint )
369         DO j=1,jjp1
370          DO i=1,iim
371           champint(i,j)=EXP(champint(i,j))
372          ENDDO
373         ENDDO
374
375         DO j = 1, jjp1
376           DO i = 1, iim
377             IF(NINT(mask(i,j)).NE.1)  THEN
378               champint( i,j ) = 0.001
379             ENDIF
380           ENDDO
381         ENDDO
382      ELSE
383         CALL rugosite(imdep, jmdep, dlon, dlat, champ,
384     .             iim, jjp1, rlonv, rlatu, champint, mask)
385      ENDIF
386         DO j = 1,jjp1
387         DO i = 1, iim
388            champtime (i,j,l) = champint(i,j)
389         ENDDO
390         ENDDO
391200      CONTINUE
392c
393      DO l = 1, lmdep
394         timeyear(l) = timecoord(l)
395      ENDDO
396
397      PRINT 222, timeyear
398222   FORMAT(2x,' Time year ',10f6.1)
399c
400       
401      PRINT*, 'Interpolation temporelle dans l annee'
402
403      DO j = 1, jjp1
404      DO i = 1, iim
405          DO l = 1, lmdep
406            ax(l) = timeyear(l)
407            ay(l) = champtime (i,j,l)
408          ENDDO
409          CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
410          DO k = 1, 360
411            time = FLOAT(k-1)
412            CALL SPLINT(ax,ay,yder,lmdep,time,by)
413            champan(i,j,k) = by
414          ENDDO
415      ENDDO
416      ENDDO
417      DO k = 1, 360
418      DO j = 1, jjp1
419         champan(iip1,j,k) = champan(1,j,k)
420      ENDDO
421        IF ( k.EQ.10 )  THEN
422          DO j = 1, jjp1
423            CALL minmax( iip1,champan(1,j,10),chmin,chmax )
424            PRINT *,' Rugosite au temps 10 ', chmin,chmax,j
425          ENDDO
426        ENDIF
427      ENDDO
428c
429      DO k = 1, 360
430         CALL gr_dyn_fi(1,iip1,jjp1,klon,champan(1,1,k), phy_rug(1,k))
431      ENDDO
432c
433      ierr = NF_CLOSE(ncid)
434
435       DEALLOCATE( dlon      )
436       DEALLOCATE( dlon_ini  )
437       DEALLOCATE( dlat      )
438       DEALLOCATE( dlat_ini  )
439       DEALLOCATE( champ     )
440c
441c
442C Traitement de la glace oceanique
443c
444      PRINT*, 'Traitement de la glace oceanique'
445      ierr = NF_OPEN('AMIP.nc', NF_NOWRITE, ncid)
446      if (ierr.ne.0) then
447        print *, NF_STRERROR(ierr)
448        STOP
449      ENDIF
450
451      ierr = NF_INQ_VARID(ncid,'SEA_ICE',varid)
452      if (ierr.ne.0) then
453        print *, NF_STRERROR(ierr)
454        STOP
455      ENDIF
456      ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
457      if (ierr.ne.0) then
458        print *, NF_STRERROR(ierr)
459        STOP
460      ENDIF
461      ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
462      if (ierr.ne.0) then
463        print *, NF_STRERROR(ierr)
464        STOP
465      ENDIF
466      print*,'variable ', namedim, 'dimension ', imdep
467      ierr = NF_INQ_VARID(ncid,namedim,dimid)
468      if (ierr.ne.0) then
469        print *, NF_STRERROR(ierr)
470        STOP
471      ENDIF
472
473      ALLOCATE ( dlon_ini(imdep) )
474      ALLOCATE (     dlon(imdep) )
475
476#ifdef NC_DOUBLE
477      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
478#else
479      ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
480#endif
481      if (ierr.ne.0) then
482        print *, NF_STRERROR(ierr)
483        STOP
484      ENDIF
485      ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
486      if (ierr.ne.0) then
487        print *, NF_STRERROR(ierr)
488        STOP
489      ENDIF
490      print*,'variable ', namedim, jmdep
491      ierr = NF_INQ_VARID(ncid,namedim,dimid)
492      if (ierr.ne.0) then
493        print *, NF_STRERROR(ierr)
494        STOP
495      ENDIF
496
497      ALLOCATE ( dlat_ini(jmdep) )
498      ALLOCATE (     dlat(jmdep) )
499
500#ifdef NC_DOUBLE
501      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
502#else
503      ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
504#endif
505      if (ierr.ne.0) then
506        print *, NF_STRERROR(ierr)
507        STOP
508      ENDIF
509      ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
510      if (ierr.ne.0) then
511        print *, NF_STRERROR(ierr)
512        STOP
513      ENDIF
514      print*,'variable ', namedim, lmdep
515      ierr = NF_INQ_VARID(ncid,namedim,dimid)
516      if (ierr.ne.0) then
517        print *, NF_STRERROR(ierr)
518        STOP
519      ENDIF
520#ifdef NC_DOUBLE
521      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
522#else
523      ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
524#endif
525      if (ierr.ne.0) then
526        print *, NF_STRERROR(ierr)
527        STOP
528      ENDIF
529c
530      ALLOCATE ( champ(imdep*jmdep) )
531
532      DO l = 1, lmdep
533         dimfirst(1) = 1
534         dimfirst(2) = 1
535         dimfirst(3) = l
536c
537         dimlast(1) = imdep
538         dimlast(2) = jmdep
539         dimlast(3) = 1
540c
541         PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
542#ifdef NC_DOUBLE
543         ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
544#else
545         ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
546#endif
547         if (ierr.ne.0) then
548           print *, NF_STRERROR(ierr)
549           STOP
550         ENDIF
551 
552         title = 'Sea-ice Amip '
553c
554         CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
555     .                        dlon, dlat, champ, interbar          )
556c
557      IF( oldice )  THEN
558                 CALL sea_ice(imdep, jmdep, dlon, dlat, champ,
559     .             iim, jjp1, rlonv, rlatu, champint )
560      ELSEIF ( interbar )  THEN
561       IF( l.EQ.1 )  THEN
562        WRITE(6,*) '-------------------------------------------------',
563     ,'------------------------'
564        WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
565     , ' pour Sea-ice Amip  $$$ '
566        WRITE(6,*) '-------------------------------------------------',
567     ,'------------------------'
568       ENDIF
569
570         CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
571     ,     champ, iim, jjm, rlonu, rlatv, jjp1, champint )
572      ELSE
573         CALL sea_ice(imdep, jmdep, dlon, dlat, champ,
574     .             iim, jjp1, rlonv, rlatu, champint )
575      ENDIF
576         DO j = 1,jjp1
577         DO i = 1, iim
578            champtime (i,j,l) = champint(i,j)
579         ENDDO
580         ENDDO
581      ENDDO
582c
583      DO l = 1, lmdep
584         timeyear(l) = timecoord(l)
585      ENDDO
586      PRINT 222,  timeyear
587c
588      PRINT*, 'Interpolation temporelle'
589      DO j = 1, jjp1
590      DO i = 1, iim
591          DO l = 1, lmdep
592            ax(l) = timeyear(l)
593            ay(l) = champtime (i,j,l)
594          ENDDO
595          CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
596          DO k = 1, 360
597            time = FLOAT(k-1)
598            CALL SPLINT(ax,ay,yder,lmdep,time,by)
599            champan(i,j,k) = by
600          ENDDO
601      ENDDO
602      ENDDO
603      DO k = 1, 360
604      DO j = 1, jjp1
605         champan(iip1, j, k) = champan(1, j, k)
606      ENDDO
607        IF ( k.EQ.10 )  THEN
608          DO j = 1, jjp1
609            CALL minmax( iip1,champan(1,j,10),chmin,chmax )
610            PRINT *,' Sea ice au temps 10 ', chmin,chmax,j
611          ENDDO
612        ENDIF
613      ENDDO
614c
615
616      DO k = 1, 360
617         CALL gr_dyn_fi(1, iip1, jjp1, klon,
618     .                  champan(1,1,k), phy_ice)
619        IF ( newlmt) THEN
620
621CPB  en attendant de mettre fraction de terre
622c
623          WHERE(phy_ice(1:klon) .GE. 1.) phy_ice(1 : klon) = 1.
624          WHERE(phy_ice(1:klon) .LT. EPSFRA) phy_ice(1 : klon) = 0.
625c
626          IF (fracterre ) THEN
627c            WRITE(*,*) 'passe dans cas fracterre'
628            pctsrf_t(:,is_ter,k) = pctsrf(:,is_ter)
629            pctsrf_t(:,is_lic,k) = pctsrf(:,is_lic)
630            pctsrf_t(1:klon,is_sic,k) =   phy_ice(1:klon)
631     $            - pctsrf_t(1:klon,is_lic,k)
632c Il y a des cas ou il y a de la glace dans landiceref et pas dans AMIP
633            WHERE (pctsrf_t(1:klon,is_sic,k) .LE. 0)
634              pctsrf_t(1:klon,is_sic,k) = 0.
635            END WHERE
636            WHERE( 1. - zmasq(1:klon) .LT. EPSFRA)
637              pctsrf_t(1:klon,is_sic,k) = 0.
638              pctsrf_t(1:klon,is_oce,k) = 0.
639            END WHERE
640            DO i = 1, klon
641              IF ( 1. - zmasq(i) .GT. EPSFRA) THEN
642                IF ( pctsrf_t(i,is_sic,k) .GE. 1 - zmasq(i)) THEN
643                  pctsrf_t(i,is_sic,k) = 1 - zmasq(i)
644                  pctsrf_t(i,is_oce,k) = 0.
645                ELSE
646                  pctsrf_t(i,is_oce,k) = 1 - zmasq(i)
647     $                    - pctsrf_t(i,is_sic,k)
648                  IF (pctsrf_t(i,is_oce,k) .LT. EPSFRA) THEN
649                    pctsrf_t(i,is_oce,k) = 0.
650                    pctsrf_t(i,is_sic,k) = 1 - zmasq(i)
651                  ENDIF
652                ENDIF
653              ENDIF 
654              if (pctsrf_t(i,is_oce,k) .lt. 0.) then
655                WRITE(*,*) 'pb sous maille au point : i,k '
656     $              , i,k,pctsrf_t(:,is_oce,k)
657              ENDIF
658              IF ( abs( pctsrf_t(i, is_ter,k) + pctsrf_t(i, is_lic,k) +
659     $          pctsrf_t(i, is_oce,k) + pctsrf_t(i, is_sic,k)  - 1.)
660     $            .GT. EPSFRA) THEN
661                  WRITE(*,*) 'physiq : pb sous surface au point ', i,
662     $                pctsrf_t(i, 1 : nbsrf,k), phy_ice(i)
663              ENDIF
664            END DO
665          ELSE
666            DO i = 1, klon
667              pctsrf_t(i,is_ter,k) = pctsrf(i,is_ter)
668              IF (NINT(pctsrf(i,is_ter)).EQ.1 ) THEN
669                pctsrf_t(i,is_sic,k) = 0.
670                pctsrf_t(i,is_oce,k) = 0.                 
671                IF(phy_ice(i) .GE. 1.e-5) THEN
672                  pctsrf_t(i,is_lic,k) = phy_ice(i)
673                  pctsrf_t(i,is_ter,k) = pctsrf_t(i,is_ter,k)
674     .                                   - pctsrf_t(i,is_lic,k)
675                ELSE
676                  pctsrf_t(i,is_lic,k) = 0.
677                ENDIF
678              ELSE
679                pctsrf_t(i,is_lic,k) = 0.
680                IF(phy_ice(i) .GE. 1.e-5) THEN
681                  pctsrf_t(i,is_ter,k) = 0.
682                  pctsrf_t(i,is_sic,k) = phy_ice(i)
683                  pctsrf_t(i,is_oce,k) = 1. - pctsrf_t(i,is_sic,k)
684                ELSE
685                  pctsrf_t(i,is_sic,k) = 0.
686                  pctsrf_t(i,is_oce,k) = 1.                     
687                ENDIF
688              ENDIF
689              verif = pctsrf_t(i,is_ter,k) +
690     .                pctsrf_t(i,is_oce,k) +
691     .                pctsrf_t(i,is_sic,k) +
692     .                pctsrf_t(i,is_lic,k)
693              IF ( verif .LT. 1. - 1.e-5 .OR.
694     $             verif .GT. 1 + 1.e-5) THEN 
695                WRITE(*,*) 'pb sous maille au point : i,k,verif '
696     $                    , i,k,verif
697              ENDIF
698            END DO
699          ENDIF
700        ELSE 
701          DO i = 1, klon
702            phy_nat(i,k) = phy_nat0(i)
703            IF ( (phy_ice(i) - 0.5).GE.1.e-5 ) THEN
704              IF (NINT(phy_nat0(i)).EQ.0) THEN
705                phy_nat(i,k) = 3.0
706              ELSE
707                phy_nat(i,k) = 2.0
708              ENDIF
709            ENDIF
710            IF( NINT(phy_nat(i,k)).EQ.0 ) THEN
711              IF ( phy_rug(i,k).NE.0.001 ) phy_rug(i,k) = 0.001
712            ENDIF
713          END DO
714        ENDIF
715      ENDDO
716c
717
718      ierr = NF_CLOSE(ncid)
719c
720       DEALLOCATE( dlon      )
721       DEALLOCATE( dlon_ini  )
722       DEALLOCATE( dlat      )
723       DEALLOCATE( dlat_ini  )
724       DEALLOCATE( champ     )
725
726477    continue
727c
728C Traitement de la sst
729c
730      PRINT*, 'Traitement de la sst'
731      ierr = NF_OPEN('AMIP.nc', NF_NOWRITE, ncid)
732      if (ierr.ne.0) then
733        print *, NF_STRERROR(ierr)
734        STOP
735      ENDIF
736
737      ierr = NF_INQ_VARID(ncid,'SST',varid)
738      if (ierr.ne.0) then
739        print *, NF_STRERROR(ierr)
740        STOP
741      ENDIF
742      ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
743      if (ierr.ne.0) then
744        print *, NF_STRERROR(ierr)
745        STOP
746      ENDIF
747      ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
748      if (ierr.ne.0) then
749        print *, NF_STRERROR(ierr)
750        STOP
751      ENDIF
752      print*,'variable SST ', namedim,'dimension ', imdep
753      ierr = NF_INQ_VARID(ncid,namedim,dimid)
754      if (ierr.ne.0) then
755        print *, NF_STRERROR(ierr)
756        STOP
757      ENDIF
758 
759      ALLOCATE( dlon_ini(imdep) )
760      ALLOCATE(     dlon(imdep) )
761
762#ifdef NC_DOUBLE
763      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
764#else
765      ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
766#endif
767
768      if (ierr.ne.0) then
769        print *, NF_STRERROR(ierr)
770        STOP
771      ENDIF
772      ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
773      if (ierr.ne.0) then
774        print *, NF_STRERROR(ierr)
775        STOP
776      ENDIF
777      print*,'variable SST ', namedim, 'dimension ', jmdep
778      ierr = NF_INQ_VARID(ncid,namedim,dimid)
779      if (ierr.ne.0) then
780        print *, NF_STRERROR(ierr)
781        STOP
782      ENDIF
783
784      ALLOCATE( dlat_ini(jmdep) )
785      ALLOCATE(     dlat(jmdep) )
786
787#ifdef NC_DOUBLE
788      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
789#else
790      ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
791#endif
792      if (ierr.ne.0) then
793        print *, NF_STRERROR(ierr)
794        STOP
795      ENDIF
796      ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
797      if (ierr.ne.0) then
798        print *, NF_STRERROR(ierr)
799        STOP
800      ENDIF
801      print*,'variable ', namedim, 'dimension ', lmdep
802      ierr = NF_INQ_VARID(ncid,namedim,dimid)
803      if (ierr.ne.0) then
804        print *, NF_STRERROR(ierr)
805        STOP
806      ENDIF
807#ifdef NC_DOUBLE
808      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
809#else
810      ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
811#endif
812      if (ierr.ne.0) then
813        print *, NF_STRERROR(ierr)
814        STOP
815      ENDIF
816
817       ALLOCATE( champ(imdep*jmdep) )
818       IF( extrap )   THEN
819         ALLOCATE ( work(imdep,jmdep) )
820       ENDIF
821c
822      DO l = 1, lmdep
823         dimfirst(1) = 1
824         dimfirst(2) = 1
825         dimfirst(3) = l
826c
827         dimlast(1) = imdep
828         dimlast(2) = jmdep
829         dimlast(3) = 1
830c
831         PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
832#ifdef NC_DOUBLE
833         ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
834#else
835         ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
836#endif
837         if (ierr.ne.0) then
838           print *, NF_STRERROR(ierr)
839           STOP
840         ENDIF
841
842         title='Sst Amip'
843c
844         CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
845     .                            dlon, dlat, champ, interbar     )
846       IF ( extrap )  THEN
847        CALL extrapol(champ, imdep, jmdep, 999999.,.TRUE.,.TRUE.,2,work)
848       ENDIF
849c
850
851      IF ( interbar )  THEN
852        IF( l.EQ.1 )  THEN
853         WRITE(6,*) '-------------------------------------------------',
854     ,'------------------------'
855         WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
856     , ' pour la Sst Amip $$$ '
857         WRITE(6,*) '-------------------------------------------------',
858     ,'------------------------'
859        ENDIF
860       CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
861     , champ, iim, jjm, rlonu, rlatv, jjp1, champint )
862      ELSE
863       CALL grille_m(imdep, jmdep, dlon, dlat, champ,
864     .          iim, jjp1, rlonv, rlatu, champint   )
865      ENDIF
866
867         DO j = 1,jjp1
868         DO i = 1, iim
869            champtime (i,j,l) = champint(i,j)
870         ENDDO
871         ENDDO
872      ENDDO
873c
874      DO l = 1, lmdep
875         timeyear(l) = timecoord(l)
876      ENDDO
877      print 222,  timeyear
878c
879C interpolation temporelle
880      DO j = 1, jjp1
881      DO i = 1, iim
882          DO l = 1, lmdep
883            ax(l) = timeyear(l)
884            ay(l) = champtime (i,j,l)
885          ENDDO
886          CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
887          DO k = 1, 360
888            time = FLOAT(k-1)
889            CALL SPLINT(ax,ay,yder,lmdep,time,by)
890            champan(i,j,k) = by
891          ENDDO
892      ENDDO
893      ENDDO
894      DO k = 1, 360
895      DO j = 1, jjp1
896         champan(iip1,j,k) = champan(1,j,k)
897      ENDDO
898        IF ( k.EQ.10 )  THEN
899          DO j = 1, jjp1
900            CALL minmax( iip1,champan(1,j,10),chmin,chmax )
901            PRINT *,' SST au temps 10 ', chmin,chmax,j
902          ENDDO
903        ENDIF
904      ENDDO
905c
906      DO k = 1, 360
907         CALL gr_dyn_fi(1, iip1, jjp1, klon,
908     .                  champan(1,1,k), phy_sst(1,k))
909      ENDDO
910c
911      ierr = NF_CLOSE(ncid)
912c
913       DEALLOCATE( dlon      )
914       DEALLOCATE( dlon_ini  )
915       DEALLOCATE( dlat      )
916       DEALLOCATE( dlat_ini  )
917       DEALLOCATE( champ     )
918c
919C Traitement de l'albedo
920c
921      PRINT*, 'Traitement de l albedo'
922      ierr = NF_OPEN('Albedo.nc', NF_NOWRITE, ncid)
923      if (ierr.ne.0) then
924        print *, NF_STRERROR(ierr)
925        STOP
926      ENDIF
927      ierr = NF_INQ_VARID(ncid,'ALBEDO',varid)
928      if (ierr.ne.0) then
929        print *, NF_STRERROR(ierr)
930        STOP
931      ENDIF
932      ierr = NF_INQ_VARDIMID (ncid,varid,ndimid)
933      if (ierr.ne.0) then
934        print *, NF_STRERROR(ierr)
935        STOP
936      ENDIF
937      ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep)
938      if (ierr.ne.0) then
939        print *, NF_STRERROR(ierr)
940        STOP
941      ENDIF
942      print*,'variable ', namedim, 'dimension ', imdep
943      ierr = NF_INQ_VARID(ncid,namedim,dimid)
944      if (ierr.ne.0) then
945        print *, NF_STRERROR(ierr)
946        STOP
947      ENDIF
948
949      ALLOCATE ( dlon_ini(imdep) )
950      ALLOCATE (     dlon(imdep) )
951
952#ifdef NC_DOUBLE
953      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlon_ini)
954#else
955      ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_ini)
956#endif
957      if (ierr.ne.0) then
958        print *, NF_STRERROR(ierr)
959        STOP
960      ENDIF
961      ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep)
962      if (ierr.ne.0) then
963        print *, NF_STRERROR(ierr)
964        STOP
965      ENDIF
966      print*,'variable ', namedim, 'dimension ', jmdep
967      ierr = NF_INQ_VARID(ncid,namedim,dimid)
968      if (ierr.ne.0) then
969        print *, NF_STRERROR(ierr)
970        STOP
971      ENDIF
972
973      ALLOCATE ( dlat_ini(jmdep) )
974      ALLOCATE (     dlat(jmdep) )
975
976#ifdef NC_DOUBLE
977      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,dlat_ini)
978#else
979      ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_ini)
980#endif
981      if (ierr.ne.0) then
982        print *, NF_STRERROR(ierr)
983        STOP
984      ENDIF
985      ierr = NF_INQ_DIM(ncid,ndimid(3), namedim, lmdep)
986      if (ierr.ne.0) then
987        print *, NF_STRERROR(ierr)
988        STOP
989      ENDIF
990      print*,'variable ', namedim, 'dimension ', lmdep
991      ierr = NF_INQ_VARID(ncid,namedim,dimid)
992      if (ierr.ne.0) then
993        print *, NF_STRERROR(ierr)
994        STOP
995      ENDIF
996#ifdef NC_DOUBLE
997      ierr = NF_GET_VAR_DOUBLE(ncid,dimid,timecoord)
998#else
999      ierr = NF_GET_VAR_REAL(ncid,dimid,timecoord)
1000#endif
1001      if (ierr.ne.0) then
1002        print *, NF_STRERROR(ierr)
1003        STOP
1004      ENDIF
1005c
1006      ALLOCATE ( champ(imdep*jmdep) )
1007
1008      DO l = 1, lmdep
1009         dimfirst(1) = 1
1010         dimfirst(2) = 1
1011         dimfirst(3) = l
1012c
1013         dimlast(1) = imdep
1014         dimlast(2) = jmdep
1015         dimlast(3) = 1
1016c
1017         PRINT*,'Lecture temporelle et int. horizontale ',l,timecoord(l)
1018#ifdef NC_DOUBLE
1019         ierr = NF_GET_VARA_DOUBLE(ncid,varid,dimfirst,dimlast,champ)
1020#else
1021         ierr = NF_GET_VARA_REAL(ncid,varid,dimfirst,dimlast,champ)
1022#endif
1023         if (ierr.ne.0) then
1024           print *, NF_STRERROR(ierr)
1025           STOP
1026         ENDIF
1027
1028         title='Albedo Amip'
1029c
1030         CALL conf_dat2d(title,imdep, jmdep, dlon_ini, dlat_ini,
1031     .                            dlon, dlat, champ, interbar      )
1032c
1033c
1034      IF ( interbar )  THEN
1035        IF( l.EQ.1 )  THEN
1036         WRITE(6,*) '-------------------------------------------------',
1037     ,'------------------------'
1038         WRITE(6,*) '$$$ Utilisation de l interpolation barycentrique ',
1039     , ' pour l Albedo Amip $$$ '
1040         WRITE(6,*) '-------------------------------------------------',
1041     ,'------------------------'
1042        ENDIF
1043
1044       CALL inter_barxy ( imdep,jmdep -1,dlon, dlat ,
1045     , champ, iim, jjm, rlonu, rlatv, jjp1, champint )
1046      ELSE
1047       CALL grille_m(imdep, jmdep, dlon, dlat, champ,
1048     .          iim, jjp1, rlonv, rlatu, champint   )
1049      ENDIF
1050c
1051         DO j = 1,jjp1
1052         DO i = 1, iim
1053            champtime (i, j, l) = champint(i, j)
1054         ENDDO
1055         ENDDO
1056      ENDDO
1057c
1058      DO l = 1, lmdep
1059         timeyear(l) = timecoord(l)
1060      ENDDO
1061      print 222,  timeyear
1062c
1063C interpolation temporelle
1064      DO j = 1, jjp1
1065      DO i = 1, iim
1066          DO l = 1, lmdep
1067            ax(l) = timeyear(l)
1068            ay(l) = champtime (i, j, l)
1069          ENDDO
1070          CALL SPLINE(ax,ay,lmdep,1.e30,1.e30,yder)
1071          DO k = 1, 360
1072            time = FLOAT(k-1)
1073            CALL SPLINT(ax,ay,yder,lmdep,time,by)
1074            champan(i,j,k) = by
1075          ENDDO
1076      ENDDO
1077      ENDDO
1078      DO k = 1, 360
1079      DO j = 1, jjp1
1080         champan(iip1, j, k) = champan(1, j, k)
1081      ENDDO
1082        IF ( k.EQ.10 )  THEN
1083          DO j = 1, jjp1
1084            CALL minmax( iip1,champan(1,j,10),chmin,chmax )
1085            PRINT *,' Albedo au temps 10 ', chmin,chmax,j
1086          ENDDO
1087        ENDIF
1088      ENDDO
1089c
1090      DO k = 1, 360
1091         CALL gr_dyn_fi(1, iip1, jjp1, klon,
1092     .                  champan(1,1,k), phy_alb(1,k))
1093      ENDDO
1094c
1095      ierr = NF_CLOSE(ncid)
1096c
1097c
1098      DO k = 1, 360
1099      DO i = 1, klon
1100         phy_bil(i,k) = 0.0
1101      ENDDO
1102      ENDDO
1103c
1104      PRINT*, 'Ecriture du fichier limit'
1105c
1106      ierr = NF_CREATE ("limit.nc", NF_CLOBBER, nid)
1107c
1108      ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 30,
1109     .                       "Fichier conditions aux limites")
1110      ierr = NF_DEF_DIM (nid, "points_physiques", klon, ndim)
1111      ierr = NF_DEF_DIM (nid, "time", NF_UNLIMITED, ntim)
1112c
1113      dims(1) = ndim
1114      dims(2) = ntim
1115c
1116      ierr = NF_DEF_VAR (nid, "TEMPS", NF_FLOAT, 1,ntim, id_tim)
1117      ierr = NF_PUT_ATT_TEXT (nid, id_tim, "title", 17,
1118     .                        "Jour dans l annee")
1119      IF (newlmt) THEN
1120c
1121        ierr = NF_DEF_VAR (nid, "FOCE", NF_FLOAT, 2,dims, id_FOCE)
1122        ierr = NF_PUT_ATT_TEXT (nid, id_FOCE, "title", 14,
1123     .                      "Fraction ocean")
1124c
1125        ierr = NF_DEF_VAR (nid, "FSIC", NF_FLOAT, 2,dims, id_FSIC)
1126        ierr = NF_PUT_ATT_TEXT (nid, id_FSIC, "title", 21,
1127     .                      "Fraction glace de mer")
1128c
1129        ierr = NF_DEF_VAR (nid, "FTER", NF_FLOAT, 2,dims, id_FTER)
1130        ierr = NF_PUT_ATT_TEXT (nid, id_FTER, "title", 14,
1131     .                      "Fraction terre")
1132c
1133        ierr = NF_DEF_VAR (nid, "FLIC", NF_FLOAT, 2,dims, id_FLIC)
1134        ierr = NF_PUT_ATT_TEXT (nid, id_FLIC, "title", 17,
1135     .                      "Fraction land ice")
1136c
1137      ELSE
1138        ierr = NF_DEF_VAR (nid, "NAT", NF_FLOAT, 2,dims, id_NAT)
1139        ierr = NF_PUT_ATT_TEXT (nid, id_NAT, "title", 23,
1140     .                      "Nature du sol (0,1,2,3)")
1141      ENDIF
1142      ierr = NF_DEF_VAR (nid, "SST", NF_FLOAT, 2,dims, id_SST)
1143      ierr = NF_PUT_ATT_TEXT (nid, id_SST, "title", 35,
1144     .                      "Temperature superficielle de la mer")
1145      ierr = NF_DEF_VAR (nid, "BILS", NF_FLOAT, 2,dims, id_BILS)
1146      ierr = NF_PUT_ATT_TEXT (nid, id_BILS, "title", 32,
1147     .                        "Reference flux de chaleur au sol")
1148      ierr = NF_DEF_VAR (nid, "ALB", NF_FLOAT, 2,dims, id_ALB)
1149      ierr = NF_PUT_ATT_TEXT (nid, id_ALB, "title", 19,
1150     .                        "Albedo a la surface")
1151      ierr = NF_DEF_VAR (nid, "RUG", NF_FLOAT, 2,dims, id_RUG)
1152      ierr = NF_PUT_ATT_TEXT (nid, id_RUG, "title", 8,
1153     .                        "Rugosite")
1154c
1155      ierr = NF_ENDDEF(nid)
1156c
1157      DO k = 1, 360
1158c
1159      debut(1) = 1
1160      debut(2) = k
1161      epais(1) = klon
1162      epais(2) = 1
1163c
1164#ifdef NC_DOUBLE
1165      ierr = NF_PUT_VAR1_DOUBLE (nid,id_tim,k,DBLE(k))
1166c
1167      IF (newlmt ) THEN
1168          ierr = NF_PUT_VARA_DOUBLE (nid,id_FOCE,debut,epais
1169     $        ,pctsrf_t(1,is_oce,k))
1170          ierr = NF_PUT_VARA_DOUBLE (nid,id_FSIC,debut,epais
1171     $        ,pctsrf_t(1,is_sic,k))
1172          ierr = NF_PUT_VARA_DOUBLE (nid,id_FTER,debut,epais
1173     $        ,pctsrf_t(1,is_ter,k))
1174          ierr = NF_PUT_VARA_DOUBLE (nid,id_FLIC,debut,epais
1175     $        ,pctsrf_t(1,is_lic,k))
1176      ELSE
1177          ierr = NF_PUT_VARA_DOUBLE (nid,id_NAT,debut,epais
1178     $        ,phy_nat(1,k))
1179      ENDIF
1180c
1181      ierr = NF_PUT_VARA_DOUBLE (nid,id_SST,debut,epais,phy_sst(1,k))
1182      ierr = NF_PUT_VARA_DOUBLE (nid,id_BILS,debut,epais,phy_bil(1,k))
1183      ierr = NF_PUT_VARA_DOUBLE (nid,id_ALB,debut,epais,phy_alb(1,k))
1184      ierr = NF_PUT_VARA_DOUBLE (nid,id_RUG,debut,epais,phy_rug(1,k))
1185#else
1186      ierr = NF_PUT_VAR1_REAL (nid,id_tim,k,FLOAT(k))
1187      IF (newlmt ) THEN
1188          ierr = NF_PUT_VARA_REAL (nid,id_FOCE,debut,epais
1189     $        ,pctsrf_t(1,is_oce,k))
1190          ierr = NF_PUT_VARA_REAL (nid,id_FSIC,debut,epais
1191     $        ,pctsrf_t(1,is_sic,k))
1192          ierr = NF_PUT_VARA_REAL (nid,id_FTER,debut,epais
1193     $        ,pctsrf_t(1,is_ter,k))
1194          ierr = NF_PUT_VARA_REAL (nid,id_FLIC,debut,epais
1195     $        ,pctsrf_t(1,is_lic,k))
1196      ELSE
1197          ierr = NF_PUT_VARA_REAL (nid,id_NAT,debut,epais
1198     $        ,phy_nat(1,k))
1199      ENDIF
1200      ierr = NF_PUT_VARA_REAL (nid,id_SST,debut,epais,phy_sst(1,k))
1201      ierr = NF_PUT_VARA_REAL (nid,id_BILS,debut,epais,phy_bil(1,k))
1202      ierr = NF_PUT_VARA_REAL (nid,id_ALB,debut,epais,phy_alb(1,k))
1203      ierr = NF_PUT_VARA_REAL (nid,id_RUG,debut,epais,phy_rug(1,k))
1204#endif
1205c
1206      ENDDO
1207c
1208      ierr = NF_CLOSE(nid)
1209c
1210      STOP
1211      END
Note: See TracBrowser for help on using the repository browser.