My Project
 All Classes Files Functions Variables Macros
dynredem_p.F
Go to the documentation of this file.
1 !
2 ! $Id: dynredem_p.F 1635 2012-07-12 11:37:16Z lguez $
3 !
4 c
5  SUBROUTINE dynredem0_p(fichnom,iday_end,phis)
6 #ifdef CPP_IOIPSL
7  USE ioipsl
8 #endif
9  USE parallel
10  USE infotrac
11  use netcdf95, only: NF95_PUT_VAR
12 
13  IMPLICIT NONE
14 c=======================================================================
15 c Ecriture du fichier de redemarrage sous format NetCDF (initialisation)
16 c=======================================================================
17 c Declarations:
18 c -------------
19 #include "dimensions.h"
20 #include "paramet.h"
21 #include "comconst.h"
22 #include "comvert.h"
23 #include "comgeom2.h"
24 #include "temps.h"
25 #include "ener.h"
26 #include "logic.h"
27 #include "netcdf.inc"
28 #include "description.h"
29 #include "serre.h"
30 
31 c Arguments:
32 c ----------
33  INTEGER iday_end
34  REAL phis(iip1, jjp1)
35  CHARACTER*(*) fichnom
36 
37 c Local:
38 c ------
39  INTEGER iq,l
40  INTEGER length
41  parameter(length = 100)
42  REAL tab_cntrl(length) ! tableau des parametres du run
43  INTEGER ierr
44  character*20 modname
45  character*80 abort_message
46 
47 c Variables locales pour NetCDF:
48 c
49  INTEGER dims2(2), dims3(3), dims4(4)
50  INTEGER idim_index
51  INTEGER idim_rlonu, idim_rlonv, idim_rlatu, idim_rlatv
52  INTEGER idim_s, idim_sig
53  INTEGER idim_tim
54  INTEGER nid,nvarid
55 
56  REAL zan0,zjulian,hours
57  INTEGER yyears0,jjour0, mmois0
58  character*30 unites
59 
60 
61 c-----------------------------------------------------------------------
62  if (mpi_rank==0) then
63 
64  modname='dynredem0_p'
65 
66 #ifdef CPP_IOIPSL
67  call ymds2ju(annee_ref, 1, iday_end, 0.0, zjulian)
68  call ju2ymds(zjulian, yyears0, mmois0, jjour0, hours)
69 #else
70 ! set yyears0, mmois0, jjour0 to 0,1,1 (hours is not used)
71  yyears0=0
72  mmois0=1
73  jjour0=1
74 #endif
75 
76  DO l=1,length
77  tab_cntrl(l) = 0.
78  ENDDO
79  tab_cntrl(1) = REAL(iim)
80  tab_cntrl(2) = REAL(jjm)
81  tab_cntrl(3) = REAL(llm)
82  tab_cntrl(4) = REAL(day_ref)
83  tab_cntrl(5) = REAL(annee_ref)
84  tab_cntrl(6) = rad
85  tab_cntrl(7) = omeg
86  tab_cntrl(8) = g
87  tab_cntrl(9) = cpp
88  tab_cntrl(10) = kappa
89  tab_cntrl(11) = daysec
90  tab_cntrl(12) = dtvr
91  tab_cntrl(13) = etot0
92  tab_cntrl(14) = ptot0
93  tab_cntrl(15) = ztot0
94  tab_cntrl(16) = stot0
95  tab_cntrl(17) = ang0
96  tab_cntrl(18) = pa
97  tab_cntrl(19) = preff
98 c
99 c ..... parametres pour le zoom ......
100 
101  tab_cntrl(20) = clon
102  tab_cntrl(21) = clat
103  tab_cntrl(22) = grossismx
104  tab_cntrl(23) = grossismy
105 c
106  IF ( fxyhypb ) THEN
107  tab_cntrl(24) = 1.
108  tab_cntrl(25) = dzoomx
109  tab_cntrl(26) = dzoomy
110  tab_cntrl(27) = 0.
111  tab_cntrl(28) = taux
112  tab_cntrl(29) = tauy
113  ELSE
114  tab_cntrl(24) = 0.
115  tab_cntrl(25) = dzoomx
116  tab_cntrl(26) = dzoomy
117  tab_cntrl(27) = 0.
118  tab_cntrl(28) = 0.
119  tab_cntrl(29) = 0.
120  IF( ysinus ) tab_cntrl(27) = 1.
121  ENDIF
122 
123  tab_cntrl(30) = REAL(iday_end)
124  tab_cntrl(31) = REAL(itau_dyn + itaufin)
125 c start_time: start_time of simulation (not necessarily 0.)
126  tab_cntrl(32) = start_time
127 c
128 c .........................................................
129 c
130 c Creation du fichier:
131 c
132  ierr = nf_create(fichnom, nf_clobber, nid)
133  IF (ierr.NE.nf_noerr) THEN
134  WRITE(6,*)" Pb d ouverture du fichier "//fichnom
135  WRITE(6,*)' ierr = ', ierr
136  CALL abort
137  ENDIF
138 c
139 c Preciser quelques attributs globaux:
140 c
141  ierr = nf_put_att_text(nid, nf_global, "title", 27,
142  . "Fichier demmarage dynamique")
143 c
144 c Definir les dimensions du fichiers:
145 c
146  ierr = nf_def_dim(nid, "index", length, idim_index)
147  ierr = nf_def_dim(nid, "rlonu", iip1, idim_rlonu)
148  ierr = nf_def_dim(nid, "rlatu", jjp1, idim_rlatu)
149  ierr = nf_def_dim(nid, "rlonv", iip1, idim_rlonv)
150  ierr = nf_def_dim(nid, "rlatv", jjm, idim_rlatv)
151  ierr = nf_def_dim(nid, "sigs", llm, idim_s)
152  ierr = nf_def_dim(nid, "sig", llmp1, idim_sig)
153  ierr = nf_def_dim(nid, "temps", nf_unlimited, idim_tim)
154 c
155  ierr = nf_enddef(nid) ! sortir du mode de definition
156 c
157 c Definir et enregistrer certains champs invariants:
158 c
159  ierr = nf_redef(nid)
160 cIM 220306 BEG
161 #ifdef NC_DOUBLE
162  ierr = nf_def_var(nid,"controle",nf_double,1,idim_index,nvarid)
163 #else
164  ierr = nf_def_var(nid,"controle",nf_float,1,idim_index,nvarid)
165 #endif
166 cIM 220306 END
167  ierr = nf_put_att_text(nid, nvarid, "title", 22,
168  . "Parametres de controle")
169  ierr = nf_enddef(nid)
170  call nf95_put_var(nid,nvarid,tab_cntrl)
171 c
172  ierr = nf_redef(nid)
173 cIM 220306 BEG
174 #ifdef NC_DOUBLE
175  ierr = nf_def_var(nid,"rlonu",nf_double,1,idim_rlonu,nvarid)
176 #else
177  ierr = nf_def_var(nid,"rlonu",nf_float,1,idim_rlonu,nvarid)
178 #endif
179 cIM 220306 END
180  ierr = nf_put_att_text(nid, nvarid, "title", 23,
181  . "Longitudes des points U")
182  ierr = nf_enddef(nid)
183  call nf95_put_var(nid,nvarid,rlonu)
184 c
185  ierr = nf_redef(nid)
186 cIM 220306 BEG
187 #ifdef NC_DOUBLE
188  ierr = nf_def_var(nid,"rlatu",nf_double,1,idim_rlatu,nvarid)
189 #else
190  ierr = nf_def_var(nid,"rlatu",nf_float,1,idim_rlatu,nvarid)
191 #endif
192 cIM 220306 END
193  ierr = nf_put_att_text(nid, nvarid, "title", 22,
194  . "Latitudes des points U")
195  ierr = nf_enddef(nid)
196  call nf95_put_var(nid,nvarid,rlatu)
197 c
198  ierr = nf_redef(nid)
199 cIM 220306 BEG
200 #ifdef NC_DOUBLE
201  ierr = nf_def_var(nid,"rlonv",nf_double,1,idim_rlonv,nvarid)
202 #else
203  ierr = nf_def_var(nid,"rlonv",nf_float,1,idim_rlonv,nvarid)
204 #endif
205 cIM 220306 END
206  ierr = nf_put_att_text(nid, nvarid, "title", 23,
207  . "Longitudes des points V")
208  ierr = nf_enddef(nid)
209  call nf95_put_var(nid,nvarid,rlonv)
210 c
211  ierr = nf_redef(nid)
212 cIM 220306 BEG
213 #ifdef NC_DOUBLE
214  ierr = nf_def_var(nid,"rlatv",nf_double,1,idim_rlatv,nvarid)
215 #else
216  ierr = nf_def_var(nid,"rlatv",nf_float,1,idim_rlatv,nvarid)
217 #endif
218 cIM 220306 END
219  ierr = nf_put_att_text(nid, nvarid, "title", 22,
220  . "Latitudes des points V")
221  ierr = nf_enddef(nid)
222  call nf95_put_var(nid,nvarid,rlatv)
223 c
224  ierr = nf_redef(nid)
225 cIM 220306 BEG
226 #ifdef NC_DOUBLE
227  ierr = nf_def_var(nid,"nivsigs",nf_double,1,idim_s,nvarid)
228 #else
229  ierr = nf_def_var(nid,"nivsigs",nf_float,1,idim_s,nvarid)
230 #endif
231 cIM 220306 END
232  ierr = nf_put_att_text(nid, nvarid, "title", 28,
233  . "Numero naturel des couches s")
234  ierr = nf_enddef(nid)
235  call nf95_put_var(nid,nvarid,nivsigs)
236 c
237  ierr = nf_redef(nid)
238 cIM 220306 BEG
239 #ifdef NC_DOUBLE
240  ierr = nf_def_var(nid,"nivsig",nf_double,1,idim_sig,nvarid)
241 #else
242  ierr = nf_def_var(nid,"nivsig",nf_float,1,idim_sig,nvarid)
243 #endif
244 cIM 220306 END
245  ierr = nf_put_att_text(nid, nvarid, "title", 32,
246  . "Numero naturel des couches sigma")
247  ierr = nf_enddef(nid)
248  call nf95_put_var(nid,nvarid,nivsig)
249 c
250  ierr = nf_redef(nid)
251 cIM 220306 BEG
252 #ifdef NC_DOUBLE
253  ierr = nf_def_var(nid,"ap",nf_double,1,idim_sig,nvarid)
254 #else
255  ierr = nf_def_var(nid,"ap",nf_float,1,idim_sig,nvarid)
256 #endif
257 cIM 220306 END
258  ierr = nf_put_att_text(nid, nvarid, "title", 26,
259  . "Coefficient A pour hybride")
260  ierr = nf_enddef(nid)
261  call nf95_put_var(nid,nvarid,ap)
262 c
263  ierr = nf_redef(nid)
264 cIM 220306 BEG
265 #ifdef NC_DOUBLE
266  ierr = nf_def_var(nid,"bp",nf_double,1,idim_sig,nvarid)
267 #else
268  ierr = nf_def_var(nid,"bp",nf_float,1,idim_sig,nvarid)
269 #endif
270 cIM 220306 END
271  ierr = nf_put_att_text(nid, nvarid, "title", 26,
272  . "Coefficient B pour hybride")
273  ierr = nf_enddef(nid)
274  call nf95_put_var(nid,nvarid,bp)
275 c
276  ierr = nf_redef(nid)
277 cIM 220306 BEG
278 #ifdef NC_DOUBLE
279  ierr = nf_def_var(nid,"presnivs",nf_double,1,idim_s,nvarid)
280 #else
281  ierr = nf_def_var(nid,"presnivs",nf_float,1,idim_s,nvarid)
282 #endif
283 cIM 220306 END
284  ierr = nf_enddef(nid)
285  call nf95_put_var(nid,nvarid,presnivs)
286 c
287 c Coefficients de passage cov. <-> contra. <--> naturel
288 c
289  ierr = nf_redef(nid)
290  dims2(1) = idim_rlonu
291  dims2(2) = idim_rlatu
292 cIM 220306 BEG
293 #ifdef NC_DOUBLE
294  ierr = nf_def_var(nid,"cu",nf_double,2,dims2,nvarid)
295 #else
296  ierr = nf_def_var(nid,"cu",nf_float,2,dims2,nvarid)
297 #endif
298 cIM 220306 END
299  ierr = nf_put_att_text(nid, nvarid, "title", 29,
300  . "Coefficient de passage pour U")
301  ierr = nf_enddef(nid)
302  call nf95_put_var(nid,nvarid,cu)
303 c
304  ierr = nf_redef(nid)
305  dims2(1) = idim_rlonv
306  dims2(2) = idim_rlatv
307 cIM 220306 BEG
308 #ifdef NC_DOUBLE
309  ierr = nf_def_var(nid,"cv",nf_double,2,dims2,nvarid)
310 #else
311  ierr = nf_def_var(nid,"cv",nf_float,2,dims2,nvarid)
312 #endif
313 cIM 220306 END
314  ierr = nf_put_att_text(nid, nvarid, "title", 29,
315  . "Coefficient de passage pour V")
316  ierr = nf_enddef(nid)
317  call nf95_put_var(nid,nvarid,cv)
318 c
319 c Aire de chaque maille:
320 c
321  ierr = nf_redef(nid)
322  dims2(1) = idim_rlonv
323  dims2(2) = idim_rlatu
324 cIM 220306 BEG
325 #ifdef NC_DOUBLE
326  ierr = nf_def_var(nid,"aire",nf_double,2,dims2,nvarid)
327 #else
328  ierr = nf_def_var(nid,"aire",nf_float,2,dims2,nvarid)
329 #endif
330 cIM 220306 END
331  ierr = nf_put_att_text(nid, nvarid, "title", 22,
332  . "Aires de chaque maille")
333  ierr = nf_enddef(nid)
334  call nf95_put_var(nid,nvarid,aire)
335 c
336 c Geopentiel au sol:
337 c
338  ierr = nf_redef(nid)
339  dims2(1) = idim_rlonv
340  dims2(2) = idim_rlatu
341 cIM 220306 BEG
342 #ifdef NC_DOUBLE
343  ierr = nf_def_var(nid,"phisinit",nf_double,2,dims2,nvarid)
344 #else
345  ierr = nf_def_var(nid,"phisinit",nf_float,2,dims2,nvarid)
346 #endif
347 cIM 220306 END
348  ierr = nf_put_att_text(nid, nvarid, "title", 19,
349  . "Geopotentiel au sol")
350  ierr = nf_enddef(nid)
351  call nf95_put_var(nid,nvarid,phis)
352 c
353 c Definir les variables pour pouvoir les enregistrer plus tard:
354 c
355  ierr = nf_redef(nid) ! entrer dans le mode de definition
356 c
357 cIM 220306 BEG
358 #ifdef NC_DOUBLE
359  ierr = nf_def_var(nid,"temps",nf_double,1,idim_tim,nvarid)
360 #else
361  ierr = nf_def_var(nid,"temps",nf_float,1,idim_tim,nvarid)
362 #endif
363 cIM 220306 END
364  ierr = nf_put_att_text(nid, nvarid, "title", 19,
365  . "Temps de simulation")
366  write(unites,200)yyears0,mmois0,jjour0
367 200 format('days since ',i4,'-',i2.2,'-',i2.2,' 00:00:00')
368  ierr = nf_put_att_text(nid, nvarid, "units", 30,
369  . unites)
370 
371 c
372  dims4(1) = idim_rlonu
373  dims4(2) = idim_rlatu
374  dims4(3) = idim_s
375  dims4(4) = idim_tim
376 cIM 220306 BEG
377 #ifdef NC_DOUBLE
378  ierr = nf_def_var(nid,"ucov",nf_double,4,dims4,nvarid)
379 #else
380  ierr = nf_def_var(nid,"ucov",nf_float,4,dims4,nvarid)
381 #endif
382 cIM 220306 END
383  ierr = nf_put_att_text(nid, nvarid, "title", 9,
384  . "Vitesse U")
385 c
386  dims4(1) = idim_rlonv
387  dims4(2) = idim_rlatv
388  dims4(3) = idim_s
389  dims4(4) = idim_tim
390 cIM 220306 BEG
391 #ifdef NC_DOUBLE
392  ierr = nf_def_var(nid,"vcov",nf_double,4,dims4,nvarid)
393 #else
394  ierr = nf_def_var(nid,"vcov",nf_float,4,dims4,nvarid)
395 #endif
396 cIM 220306 END
397  ierr = nf_put_att_text(nid, nvarid, "title", 9,
398  . "Vitesse V")
399 c
400  dims4(1) = idim_rlonv
401  dims4(2) = idim_rlatu
402  dims4(3) = idim_s
403  dims4(4) = idim_tim
404 cIM 220306 BEG
405 #ifdef NC_DOUBLE
406  ierr = nf_def_var(nid,"teta",nf_double,4,dims4,nvarid)
407 #else
408  ierr = nf_def_var(nid,"teta",nf_float,4,dims4,nvarid)
409 #endif
410 cIM 220306 END
411  ierr = nf_put_att_text(nid, nvarid, "title", 11,
412  . "Temperature")
413 c
414  dims4(1) = idim_rlonv
415  dims4(2) = idim_rlatu
416  dims4(3) = idim_s
417  dims4(4) = idim_tim
418 
419  DO iq=1,nqtot
420 cIM 220306 BEG
421 #ifdef NC_DOUBLE
422  ierr = nf_def_var(nid,tname(iq),nf_double,4,dims4,nvarid)
423 #else
424  ierr = nf_def_var(nid,tname(iq),nf_float,4,dims4,nvarid)
425 #endif
426 cIM 220306 END
427  ierr = nf_put_att_text(nid, nvarid, "title", 12,ttext(iq))
428  ENDDO
429 c
430  dims4(1) = idim_rlonv
431  dims4(2) = idim_rlatu
432  dims4(3) = idim_s
433  dims4(4) = idim_tim
434 cIM 220306 BEG
435 #ifdef NC_DOUBLE
436  ierr = nf_def_var(nid,"masse",nf_double,4,dims4,nvarid)
437 #else
438  ierr = nf_def_var(nid,"masse",nf_float,4,dims4,nvarid)
439 #endif
440 cIM 220306 END
441  ierr = nf_put_att_text(nid, nvarid, "title", 12,
442  . "C est quoi ?")
443 c
444  dims3(1) = idim_rlonv
445  dims3(2) = idim_rlatu
446  dims3(3) = idim_tim
447 cIM 220306 BEG
448 #ifdef NC_DOUBLE
449  ierr = nf_def_var(nid,"ps",nf_double,3,dims3,nvarid)
450 #else
451  ierr = nf_def_var(nid,"ps",nf_float,3,dims3,nvarid)
452 #endif
453 cIM 220306 END
454  ierr = nf_put_att_text(nid, nvarid, "title", 15,
455  . "Pression au sol")
456 c
457  ierr = nf_enddef(nid) ! sortir du mode de definition
458  ierr = nf_close(nid) ! fermer le fichier
459 
460  print*,'iim,jjm,llm,iday_end',iim,jjm,llm,iday_end
461  print*,'rad,omeg,g,cpp,kappa',
462  , rad,omeg,g,cpp,kappa
463 
464  endif ! mpi_rank==0
465  RETURN
466  END
467  SUBROUTINE dynredem1_p(fichnom,time,
468  . vcov,ucov,teta,q,masse,ps)
469  USE parallel
470  USE infotrac
471  USE control_mod
472  use netcdf, only: NF90_get_VAR
473  use netcdf95, only: NF95_PUT_VAR
474 
475  IMPLICIT NONE
476 c=================================================================
477 c Ecriture du fichier de redemarrage sous format NetCDF
478 c=================================================================
479 #include "dimensions.h"
480 #include "paramet.h"
481 #include "description.h"
482 #include "netcdf.inc"
483 #include "comvert.h"
484 #include "comgeom.h"
485 #include "temps.h"
486 
487 
488  INTEGER l
489  REAL vcov(iip1,jjm,llm),ucov(iip1, jjp1,llm)
490  REAL teta(iip1, jjp1,llm)
491  REAL ps(iip1, jjp1),masse(iip1, jjp1,llm)
492  REAL q(iip1, jjp1, llm, nqtot)
493  CHARACTER*(*) fichnom
494 
495  REAL time
496  INTEGER nid, nvarid, nid_trac, nvarid_trac
497  REAL trac_tmp(ip1jmp1,llm)
498  INTEGER ierr, ierr_file
499  INTEGER iq
500  INTEGER length
501  parameter(length = 100)
502  REAL tab_cntrl(length) ! tableau des parametres du run
503  character*20 modname
504  character*80 abort_message
505 c
506  INTEGER nb
507  SAVE nb
508  DATA nb / 0 /
509 
510  logical exist_file
511 
512  call gather_field(ucov,ip1jmp1,llm,0)
513  call gather_field(vcov,ip1jm,llm,0)
514  call gather_field(teta,ip1jmp1,llm,0)
515  call gather_field(masse,ip1jmp1,llm,0)
516  call gather_field(ps,ip1jmp1,1,0)
517 
518  do iq=1,nqtot
519  call gather_field(q(:,:,:,iq),ip1jmp1,llm,0)
520  enddo
521 
522 
523  if (mpi_rank==0) then
524 
525  modname = 'dynredem1'
526  ierr = nf_open(fichnom, nf_write, nid)
527  IF (ierr .NE. nf_noerr) THEN
528  print*, "Pb. d ouverture "//fichnom
529  CALL abort
530  ENDIF
531 
532 c Ecriture/extension de la coordonnee temps
533 
534  nb = nb + 1
535  ierr = nf_inq_varid(nid, "temps", nvarid)
536  IF (ierr .NE. nf_noerr) THEN
537  print *, nf_strerror(ierr)
538  abort_message='Variable temps n est pas definie'
539  CALL abort_gcm(modname,abort_message,ierr)
540  ENDIF
541  call nf95_put_var(nid,nvarid,time,start=(/nb/))
542  print*, "Enregistrement pour ", nb, time
543 
544 c
545 c Re-ecriture du tableau de controle, itaufin n'est plus defini quand
546 c on passe dans dynredem0
547  ierr = nf_inq_varid(nid, "controle", nvarid)
548  IF (ierr .NE. nf_noerr) THEN
549  abort_message="dynredem1: Le champ <controle> est absent"
550  ierr = 1
551  CALL abort_gcm(modname,abort_message,ierr)
552  ENDIF
553  ierr = nf90_get_var(nid, nvarid, tab_cntrl)
554  tab_cntrl(31) = REAL(itau_dyn + itaufin)
555  call nf95_put_var(nid,nvarid,tab_cntrl)
556 
557 c Ecriture des champs
558 c
559  ierr = nf_inq_varid(nid, "ucov", nvarid)
560  IF (ierr .NE. nf_noerr) THEN
561  print*, "Variable ucov n est pas definie"
562  CALL abort
563  ENDIF
564  call nf95_put_var(nid,nvarid,ucov)
565 
566  ierr = nf_inq_varid(nid, "vcov", nvarid)
567  IF (ierr .NE. nf_noerr) THEN
568  print*, "Variable vcov n est pas definie"
569  CALL abort
570  ENDIF
571  call nf95_put_var(nid,nvarid,vcov)
572 
573  ierr = nf_inq_varid(nid, "teta", nvarid)
574  IF (ierr .NE. nf_noerr) THEN
575  print*, "Variable teta n est pas definie"
576  CALL abort
577  ENDIF
578  call nf95_put_var(nid,nvarid,teta)
579 
580  IF (type_trac == 'inca') THEN
581 ! Ajout Anne pour lecture valeurs traceurs dans un fichier start_trac.nc
582  inquire(file="start_trac.nc", exist=exist_file)
583  print *, "EXIST", exist_file
584  if (exist_file) then
585  ierr_file = nf_open("start_trac.nc", nf_nowrite,nid_trac)
586  IF (ierr_file .NE.nf_noerr) THEN
587  write(6,*)' Pb d''ouverture du fichier start_trac.nc'
588  write(6,*)' ierr = ', ierr_file
589  ENDIF
590  else
591  ierr_file = 2
592  endif
593  END IF
594 
595  do iq=1,nqtot
596 
597  IF (type_trac /= 'inca') THEN
598  ierr = nf_inq_varid(nid, tname(iq), nvarid)
599  IF (ierr .NE. nf_noerr) THEN
600  print*, "Variable tname(iq) n est pas definie"
601  CALL abort
602  ENDIF
603  call nf95_put_var(nid,nvarid,q(:,:,:,iq))
604  ELSE ! type_trac = inca
605 ! lecture de la valeur du traceur dans start_trac.nc
606  IF (ierr_file .ne. 2) THEN
607  ierr = nf_inq_varid(nid_trac, tname(iq), nvarid_trac)
608  IF (ierr .NE. nf_noerr) THEN
609  print*, tname(iq),"est absent de start_trac.nc"
610  ierr = nf_inq_varid(nid, tname(iq), nvarid)
611  IF (ierr .NE. nf_noerr) THEN
612  print*, "Variable ", tname(iq)," n est pas definie"
613  CALL abort
614  ENDIF
615  call nf95_put_var(nid,nvarid,q(:,:,:,iq))
616 
617  ELSE
618  print*, tname(iq), "est present dans start_trac.nc"
619  ierr = nf90_get_var(nid_trac, nvarid_trac, trac_tmp)
620  IF (ierr .NE. nf_noerr) THEN
621  print*, "Lecture echouee pour", tname(iq)
622  CALL abort
623  ENDIF
624  ierr = nf_inq_varid(nid, tname(iq), nvarid)
625  IF (ierr .NE. nf_noerr) THEN
626  print*, "Variable ", tname(iq)," n est pas definie"
627  CALL abort
628  ENDIF
629  call nf95_put_var(nid, nvarid, trac_tmp)
630 
631  ENDIF ! IF (ierr .NE. NF_NOERR)
632 ! fin lecture du traceur
633  ELSE ! si il n'y a pas de fichier start_trac.nc
634 ! print *, 'il n y a pas de fichier start_trac'
635  ierr = nf_inq_varid(nid, tname(iq), nvarid)
636  IF (ierr .NE. nf_noerr) THEN
637  print*, "Variable tname(iq) n est pas definie"
638  CALL abort
639  ENDIF
640  call nf95_put_var(nid,nvarid,q(:,:,:,iq))
641  ENDIF ! (ierr_file .ne. 2)
642  END IF !type_trac
643 
644  ENDDO
645 
646 
647 
648 c
649  ierr = nf_inq_varid(nid, "masse", nvarid)
650  IF (ierr .NE. nf_noerr) THEN
651  print*, "Variable masse n est pas definie"
652  CALL abort
653  ENDIF
654  call nf95_put_var(nid,nvarid,masse)
655 c
656  ierr = nf_inq_varid(nid, "ps", nvarid)
657  IF (ierr .NE. nf_noerr) THEN
658  print*, "Variable ps n est pas definie"
659  CALL abort
660  ENDIF
661  call nf95_put_var(nid,nvarid,ps)
662 
663  ierr = nf_close(nid)
664 c
665  endif ! mpi_rank==0
666 
667  RETURN
668  END
669