Context Navigation

newstart.F @ 1356

Last change on this file since 1356 was 1356, checked in by slebonnois, 10 years ago
SL: update to newstart/start2archive tools in Venus+Titan / additional diagnostics in radiative fluxes for Titan
File size: 39.2 KB

Line
1	C======================================================================
2	PROGRAM newstart
3	c=======================================================================
4	c
5	c
6	c Auteur: S. Lebonnois,
7	c a partir des newstart/start_archive/lect_start_archive martiens
8	c
9	c Derniere modif : 02/09 (ecriture des q*)
10	c 01/12 (inclusion dans svn dyn3d)
11	c
12	c Objet: Modify the grid for the initial state (LMD GCM VENUS/TITAN)
13	c ----- (from file NetCDF start_archive.nc)
14	c
15	c
16	c=======================================================================
17
18	use IOIPSL
19	USE filtreg_mod
20	USE startvar
21	USE control_mod
22	USE infotrac
23	use cpdet_mod, only: ini_cpdet,t2tpot
24	use exner_hyb_m, only: exner_hyb
25	use exner_milieu_m, only: exner_milieu
26
27	implicit none
28
29	#include "dimensions.h"
30	#include "paramet.h"
31	#include "comconst.h"
32	#include "comdissnew.h"
33	#include "comvert.h"
34	#include "comgeom2.h"
35	#include "logic.h"
36	#include "temps.h"
37	#include "ener.h"
38	#include "description.h"
39	#include "serre.h"
40	#include "dimsoil.h"
41	#include "netcdf.inc"
42
43	c-----------------------------------------------------------------------
44	c Declarations
45	c-----------------------------------------------------------------------
46
47	c Variables pour fichier "ini"
48	c------------------------------------
49	INTEGER imold,jmold,lmold,nqold,ip1jmp1old
50	INTEGER length
51	parameter (length = 100)
52	real tab_cntrl(2*length)
53	INTEGER isoil,iq,iqmax
54	CHARACTER*2 str2
55
56	c Variable histoire
57	c------------------
58	REAL vcov(iip1,jjm,llm),ucov(iip1,jjp1,llm) ! vents covariants
59	REAL teta(iip1,jjp1,llm),ps(iip1,jjp1)
60	REAL phis(iip1,jjp1) ! geopotentiel au sol
61	REAL masse(ip1jmp1,llm) ! masse de l'atmosphere
62	REAL, ALLOCATABLE, DIMENSION(:,:,:,:):: q! champs advectes
63	REAL tab_cntrl_dyn(length) ! tableau des parametres de start
64
65	c variable physique
66	c------------------
67	integer ngridmx
68	parameter (ngridmx=(2+(jjm-1)*iim - 1/jjm))
69	REAL tab_cntrl_fi(length) ! tableau des parametres de startfi
70	real rlat(ngridmx),rlon(ngridmx)
71	REAL tsurf(ngridmx),tsoil(ngridmx,nsoilmx)
72	REAL albe(ngridmx),radsol(ngridmx),sollw(ngridmx)
73	real solsw(ngridmx),fder(ngridmx)
74	real sollwdown(ngridmx),dlw(ngridmx)
75	REAL zmea(ngridmx), zstd(ngridmx)
76	REAL zsig(ngridmx), zgam(ngridmx), zthe(ngridmx)
77	REAL zpic(ngridmx), zval(ngridmx)
78	real t_fi(ngridmx,llm)
79
80	c Variable nouvelle grille naturelle au point scalaire
81	c------------------------------------------------------
82	real us(iip1,jjp1,llm),vs(iip1,jjp1,llm)
83	REAL p3d(iip1,jjp1,llm+1) ! pression aux interfaces
84	REAL phisold_newgrid(iip1,jjp1)
85	REAL T(iip1,jjp1,llm)
86	real rlonS(iip1,jjp1),rlatS(iip1,jjp1)
87	real tsurfS(iip1,jjp1),tsoilS(iip1,jjp1,nsoilmx)
88	real albeS(ip1jmp1),radsolS(ip1jmp1),sollwS(ip1jmp1)
89	real solswS(ip1jmp1),fderS(ip1jmp1)
90	real sollwdownS(ip1jmp1),dlwS(ip1jmp1)
91	real zmeaS(ip1jmp1),zstdS(ip1jmp1),zsigS(ip1jmp1)
92	real zgamS(ip1jmp1),ztheS(ip1jmp1),zpicS(ip1jmp1)
93	real zvalS(ip1jmp1)
94
95	real ptotal
96
97	c Var intermediaires : vent naturel, mais pas coord scalaire
98	c-----------------------------------------------------------
99	real vnat(iip1,jjm,llm),unat(iip1,jjp1,llm)
100
101	REAL pks(iip1,jjp1) ! exner (f pour filtre)
102	REAL pk(iip1,jjp1,llm)
103	REAL pkf(iip1,jjp1,llm)
104	REAL alpha(iip1,jjp1,llm),beta(iip1,jjp1,llm)
105
106
107	c Variable de l'ancienne grille
108	c---------------------------------------------------------
109
110	real, dimension(:), allocatable :: rlonuold, rlatvold
111	real, dimension(:), allocatable :: rlonvold, rlatuold
112	real, dimension(:), allocatable :: nivsigsold,nivsigold
113	real, dimension(:), allocatable :: apold,bpold
114	real, dimension(:), allocatable :: presnivsold
115	real, dimension(:,:,:), allocatable :: uold,vold,told
116	real, dimension(:,:,:,:), allocatable :: qold
117	real, dimension(:,:,:), allocatable :: tsoilold
118	real, dimension(:,:), allocatable :: psold,phisold
119	real, dimension(:,:), allocatable :: tsurfold
120	real, dimension(:,:), allocatable :: albeold,radsolold
121	real, dimension(:,:), allocatable :: sollwold,solswold
122	real, dimension(:,:), allocatable :: fderold
123	real, dimension(:,:), allocatable :: dlwold,sollwdownold
124	real, dimension(:,:), allocatable :: zmeaold,zstdold,zsigold
125	real, dimension(:,:), allocatable :: zgamold,ztheold,zpicold
126	real, dimension(:,:), allocatable :: zvalold
127
128	real ptotalold
129
130	c Variable intermediaires iutilise pour l'extrapolation verticale
131	c----------------------------------------------------------------
132	real, dimension(:,:,:), allocatable :: var,varp1
133
134	c divers local
135	c-----------------
136
137	integer ierr,nid,nvarid
138	INTEGER ij, l,i,j
139	character*80 fichnom
140	integer, dimension(4) :: start,counter
141	REAL phisinverse(iip1,jjp1) ! geopotentiel au sol avant inversion
142	logical topoflag,albedoflag,razvitu,razvitv
143	real albedo
144
145	c=======================================================================
146	c INITIALISATIONS DIVERSES
147	c=======================================================================
148
149	c VENUS/TITAN
150
151	iflag_trac = 1
152	c-----------------------------------------------------------------------
153	c Initialisation des traceurs
154	c ---------------------------
155	c Choix du nombre de traceurs et du schema pour l'advection
156	c dans fichier traceur.def, par default ou via INCA
157	call infotrac_init
158
159	c Allocation de la tableau q : champs advectes
160	allocate(q(iip1,jjp1,llm,nqtot))
161
162	c-----------------------------------------------------------------------
163	c Ouverture du fichier a modifier (start_archive.nc)
164	c-----------------------------------------------------------------------
165
166	write(,) 'Creation d un etat initial a partir de'
167	write(,) './start_archive.nc'
168	write(,)
169	fichnom = 'start_archive.nc'
170	ierr = NF_OPEN (fichnom, NF_NOWRITE,nid)
171	IF (ierr.NE.NF_NOERR) THEN
172	write(6,*)' Pb d''ouverture du fichier ',fichnom
173	write(6,*)' ierr = ', ierr
174	CALL ABORT
175	ENDIF
176
177	c-----------------------------------------------------------------------
178	c Lecture du tableau des parametres du run (pour la dynamique)
179	c-----------------------------------------------------------------------
180
181	write(,) 'lecture tab_cntrl START_ARCHIVE'
182	c
183	ierr = NF_INQ_VARID (nid, "controle", nvarid)
184	#ifdef NC_DOUBLE
185	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tab_cntrl)
186	#else
187	ierr = NF_GET_VAR_REAL(nid, nvarid, tab_cntrl)
188	#endif
189	c
190	write(,) 'Impression de tab_cntrl'
191	do i=1,200
192	write(,) i,tab_cntrl(i)
193	enddo
194
195	c-----------------------------------------------------------------------
196	c Initialisation des constantes
197	c-----------------------------------------------------------------------
198
199	imold = tab_cntrl(1)
200	jmold = tab_cntrl(2)
201	lmold = tab_cntrl(3)
202	day_ref = tab_cntrl(4)
203	annee_ref = tab_cntrl(5)
204	rad = tab_cntrl(6)
205	omeg = tab_cntrl(7)
206	g = tab_cntrl(8)
207	cpp = tab_cntrl(9)
208	kappa = tab_cntrl(10)
209	daysec = tab_cntrl(11)
210	dtvr = tab_cntrl(12)
211	etot0 = tab_cntrl(13)
212	ptot0 = tab_cntrl(14)
213	ztot0 = tab_cntrl(15)
214	stot0 = tab_cntrl(16)
215	ang0 = tab_cntrl(17)
216	pa = tab_cntrl(18)
217	preff = tab_cntrl(19)
218	c
219	clon = tab_cntrl(20)
220	clat = tab_cntrl(21)
221	grossismx = tab_cntrl(22)
222	grossismy = tab_cntrl(23)
223
224	IF ( tab_cntrl(24).EQ.1. ) THEN
225	fxyhypb = . TRUE .
226	dzoomx = tab_cntrl(25)
227	dzoomy = tab_cntrl(26)
228	taux = tab_cntrl(28)
229	tauy = tab_cntrl(29)
230	ELSE
231	fxyhypb = . FALSE .
232	ysinus = . FALSE .
233	IF( tab_cntrl(27).EQ.1. ) ysinus = . TRUE.
234	ENDIF
235
236	ptotalold = tab_cntrl(2*length)
237
238	write(,) "Old dimensions:"
239	write(,) "longitude: ",imold
240	write(,) "latitude: ",jmold
241	write(,) "altitude: ",lmold
242	write(,)
243
244	ip1jmp1old = (imold+1-1/imold)*(jmold+1-1/jmold)
245
246	c dans run.def
247	CALL getin('planet_type',planet_type)
248	call ini_cpdet
249
250	c=======================================================================
251	c CHANGEMENT DE CONSTANTES CONTENUES DANS tab_cntrl
252	c=======================================================================
253	c changement de la resolution dans le fichier de controle
254	tab_cntrl(1) = REAL(iim)
255	tab_cntrl(2) = REAL(jjm)
256	tab_cntrl(3) = REAL(llm)
257
258	c--------------------------------
259	c We use a specific run.def to read new parameters that need to be changed
260	c--------------------------------
261
262	c Changement de cpp:
263	CALL getin('cpp',cpp)
264	kappa = (8.314511/43.44E-3)/cpp
265	tab_cntrl(9) = cpp
266	tab_cntrl(10) = kappa
267
268	c CHANGING THE ZOOM PARAMETERS TO CHANGE THE GRID
269	CALL getin('clon',clon)
270	CALL getin('clat',clat)
271	tab_cntrl(20) = clon
272	tab_cntrl(21) = clat
273	CALL getin('grossismx',grossismx)
274	CALL getin('grossismy',grossismy)
275	tab_cntrl(22) = grossismx
276	tab_cntrl(23) = grossismy
277	CALL getin('fxyhypb',fxyhypb)
278	IF ( fxyhypb ) THEN
279	CALL getin('dzoomx',dzoomx)
280	CALL getin('dzoomy',dzoomy)
281	tab_cntrl(25) = dzoomx
282	tab_cntrl(26) = dzoomy
283	CALL getin('taux',taux)
284	CALL getin('tauy',tauy)
285	tab_cntrl(28) = taux
286	tab_cntrl(29) = tauy
287	ELSE
288	CALL getin('ysinus',ysinus)
289	tab_cntrl(27) = ysinus
290	ENDIF
291
292	c changes must be done BEFORE these lines...
293	DO l=1,length
294	tab_cntrl_dyn(l)= tab_cntrl(l)
295	tab_cntrl_fi(l) = tab_cntrl(l+length)
296	ENDDO
297
298	c-----------------------------------------------------------------------
299	c Autres initialisations
300	c-----------------------------------------------------------------------
301
302	CALL iniconst
303	CALL inigeom
304	call inifilr
305
306	c-----------------------------------------------------------------------
307	c Allocations des anciennes variables
308	c-----------------------------------------------------------------------
309
310	allocate(rlonuold(imold+1), rlatvold(jmold))
311	allocate(rlonvold(imold+1), rlatuold(jmold+1))
312	allocate(nivsigsold(lmold+1),nivsigold(lmold))
313	allocate(apold(lmold),bpold(lmold))
314	allocate(presnivsold(lmold))
315	allocate(uold(imold+1,jmold+1,lmold))
316	allocate(vold(imold+1,jmold+1,lmold))
317	allocate(told(imold+1,jmold+1,lmold))
318	allocate(qold(imold+1,jmold+1,lmold,nqtot))
319	allocate(psold(imold+1,jmold+1))
320	allocate(phisold(imold+1,jmold+1))
321	allocate(tsurfold(imold+1,jmold+1))
322	allocate(tsoilold(imold+1,jmold+1,nsoilmx))
323	allocate(albeold(imold+1,jmold+1),radsolold(imold+1,jmold+1))
324	allocate(sollwold(imold+1,jmold+1),solswold(imold+1,jmold+1))
325	allocate(fderold(imold+1,jmold+1))
326	allocate(sollwdownold(imold+1,jmold+1),dlwold(imold+1,jmold+1))
327	allocate(zmeaold(imold+1,jmold+1),zstdold(imold+1,jmold+1))
328	allocate(zsigold(imold+1,jmold+1),zgamold(imold+1,jmold+1))
329	allocate(ztheold(imold+1,jmold+1),zpicold(imold+1,jmold+1))
330	allocate(zvalold(imold+1,jmold+1))
331
332	allocate(var (imold+1,jmold+1,llm))
333	allocate(varp1 (imold+1,jmold+1,llm+1))
334
335	print*,"Initialisations OK"
336
337	c-----------------------------------------------------------------------
338	c Lecture des longitudes et latitudes
339	c-----------------------------------------------------------------------
340	c
341	ierr = NF_INQ_VARID (nid, "rlonv", nvarid)
342	IF (ierr .NE. NF_NOERR) THEN
343	PRINT*, "new_grid: Le champ <rlonv> est absent de start.nc"
344	CALL abort
345	ENDIF
346	#ifdef NC_DOUBLE
347	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlonvold)
348	#else
349	ierr = NF_GET_VAR_REAL(nid, nvarid, rlonvold)
350	#endif
351	IF (ierr .NE. NF_NOERR) THEN
352	PRINT*, "new_grid: Lecture echouee pour <rlonv>"
353	CALL abort
354	ENDIF
355	c
356	ierr = NF_INQ_VARID (nid, "rlatu", nvarid)
357	IF (ierr .NE. NF_NOERR) THEN
358	PRINT*, "new_grid: Le champ <rlatu> est absent de start.nc"
359	CALL abort
360	ENDIF
361	#ifdef NC_DOUBLE
362	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlatuold)
363	#else
364	ierr = NF_GET_VAR_REAL(nid, nvarid, rlatuold)
365	#endif
366	IF (ierr .NE. NF_NOERR) THEN
367	PRINT*, "new_grid: Lecture echouee pour <rlatu>"
368	CALL abort
369	ENDIF
370	c
371	ierr = NF_INQ_VARID (nid, "rlonu", nvarid)
372	IF (ierr .NE. NF_NOERR) THEN
373	PRINT*, "new_grid: Le champ <rlonu> est absent de start.nc"
374	CALL abort
375	ENDIF
376	#ifdef NC_DOUBLE
377	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlonuold)
378	#else
379	ierr = NF_GET_VAR_REAL(nid, nvarid, rlonuold)
380	#endif
381	IF (ierr .NE. NF_NOERR) THEN
382	PRINT*, "new_grid: Lecture echouee pour <rlonu>"
383	CALL abort
384	ENDIF
385	c
386	ierr = NF_INQ_VARID (nid, "rlatv", nvarid)
387	IF (ierr .NE. NF_NOERR) THEN
388	PRINT*, "new_grid: Le champ <rlatv> est absent de start.nc"
389	CALL abort
390	ENDIF
391	#ifdef NC_DOUBLE
392	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rlatvold)
393	#else
394	ierr = NF_GET_VAR_REAL(nid, nvarid, rlatvold)
395	#endif
396	IF (ierr .NE. NF_NOERR) THEN
397	PRINT*, "new_grid: Lecture echouee pour <rlatv>"
398	CALL abort
399	ENDIF
400	c
401
402	print*,"Lecture lon/lat OK"
403
404	c-----------------------------------------------------------------------
405	c Lecture des niveaux verticaux
406	c-----------------------------------------------------------------------
407	c
408	ierr = NF_INQ_VARID (nid, "nivsigs", nvarid)
409	IF (ierr .NE. NF_NOERR) THEN
410	PRINT*, "dynetat0: Le champ <nivsigs> est absent"
411	CALL abort
412	ENDIF
413	#ifdef NC_DOUBLE
414	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, nivsigsold)
415	#else
416	ierr = NF_GET_VAR_REAL(nid, nvarid, nivsigsold)
417	#endif
418	IF (ierr .NE. NF_NOERR) THEN
419	PRINT*, "dynetat0: Lecture echouee pour <nivsigs>"
420	CALL abort
421	ENDIF
422
423	ierr = NF_INQ_VARID (nid, "nivsig", nvarid)
424	IF (ierr .NE. NF_NOERR) THEN
425	PRINT*, "dynetat0: Le champ <nivsig> est absent"
426	CALL abort
427	ENDIF
428	#ifdef NC_DOUBLE
429	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, nivsigold)
430	#else
431	ierr = NF_GET_VAR_REAL(nid, nvarid, nivsigold)
432	#endif
433	IF (ierr .NE. NF_NOERR) THEN
434	PRINT*, "dynetat0: Lecture echouee pour <nivsig>"
435	CALL abort
436	ENDIF
437
438	ierr = NF_INQ_VARID (nid, "ap", nvarid)
439	IF (ierr .NE. NF_NOERR) THEN
440	PRINT*, "new_grid: Le champ <ap> est absent de start.nc"
441	CALL abort
442	ELSE
443	#ifdef NC_DOUBLE
444	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, apold)
445	#else
446	ierr = NF_GET_VAR_REAL(nid, nvarid, apold)
447	#endif
448	IF (ierr .NE. NF_NOERR) THEN
449	PRINT*, "new_grid: Lecture echouee pour <ap>"
450	ENDIF
451	ENDIF
452	c
453	ierr = NF_INQ_VARID (nid, "bp", nvarid)
454	IF (ierr .NE. NF_NOERR) THEN
455	PRINT*, "new_grid: Le champ <bp> est absent de start.nc"
456	CALL abort
457	ENDIF
458	#ifdef NC_DOUBLE
459	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, bpold)
460	#else
461	ierr = NF_GET_VAR_REAL(nid, nvarid, bpold)
462	#endif
463	IF (ierr .NE. NF_NOERR) THEN
464	PRINT*, "new_grid: Lecture echouee pour <bp>"
465	CALL abort
466	END IF
467
468	ierr = NF_INQ_VARID (nid, "presnivs", nvarid)
469	IF (ierr .NE. NF_NOERR) THEN
470	PRINT*, "dynetat0: Le champ <presnivs> est absent"
471	CALL abort
472	ENDIF
473	#ifdef NC_DOUBLE
474	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, presnivsold)
475	#else
476	ierr = NF_GET_VAR_REAL(nid, nvarid, presnivsold)
477	#endif
478	IF (ierr .NE. NF_NOERR) THEN
479	PRINT*, "dynetat0: Lecture echouee pour <presnivs>"
480	CALL abort
481	ENDIF
482
483	c-----------------------------------------------------------------------
484	c Lecture geopotentiel au sol
485	c-----------------------------------------------------------------------
486	c
487	ierr = NF_INQ_VARID (nid, "phisinit", nvarid)
488	IF (ierr .NE. NF_NOERR) THEN
489	PRINT*, "new_grid: Le champ <phisinit> est absent de start.nc"
490	CALL abort
491	ENDIF
492	#ifdef NC_DOUBLE
493	ierr = NF_GET_VAR_DOUBLE(nid, nvarid, phisold)
494	#else
495	ierr = NF_GET_VAR_REAL(nid, nvarid, phisold)
496	#endif
497	IF (ierr .NE. NF_NOERR) THEN
498	PRINT*, "new_grid: Lecture echouee pour <phisinit>"
499	CALL abort
500	ENDIF
501
502	print*,"Lecture niveaux et geopot OK"
503
504	c-----------------------------------------------------------------------
505	c Lecture des champs 2D ()
506	c-----------------------------------------------------------------------
507
508	start=(/1,1,1,0/)
509	counter=(/imold+1,jmold+1,1,0/)
510
511	ierr = NF_INQ_VARID (nid, "ps", nvarid)
512	IF (ierr .NE. NF_NOERR) THEN
513	PRINT*, "lect_start_archive: Le champ <ps> est absent"
514	CALL abort
515	ENDIF
516	#ifdef NC_DOUBLE
517	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,psold)
518	#else
519	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,psold)
520	#endif
521	IF (ierr .NE. NF_NOERR) THEN
522	PRINT*, "lect_start_archive: Lecture echouee pour <ps>"
523	CALL abort
524	ENDIF
525	c
526	ierr = NF_INQ_VARID (nid, "tsurf", nvarid)
527	IF (ierr .NE. NF_NOERR) THEN
528	PRINT*, "lect_start_archive: Le champ <tsurf> est absent"
529	CALL abort
530	ENDIF
531	#ifdef NC_DOUBLE
532	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,tsurfold)
533	#else
534	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,tsurfold)
535	#endif
536	IF (ierr .NE. NF_NOERR) THEN
537	PRINT*, "lect_start_archive: Lecture echouee pour <tsurf>"
538	CALL abort
539	ENDIF
540	c
541	do isoil=1,nsoilmx
542	write(str2,'(i2.2)') isoil
543	c
544	ierr = NF_INQ_VARID (nid, "Tsoil"//str2, nvarid)
545	IF (ierr .NE. NF_NOERR) THEN
546	PRINT*, "lect_start_archive:
547	. Le champ <","Tsoil"//str2,"> est absent"
548	CALL abort
549	ENDIF
550	#ifdef NC_DOUBLE
551	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,
552	. tsoilold(1,1,isoil))
553	#else
554	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,
555	. tsoilold(1,1,isoil))
556	#endif
557	IF (ierr .NE. NF_NOERR) THEN
558	PRINT*, "lect_start_archive:
559	. Lecture echouee pour <","Tsoil"//str2,">"
560	CALL abort
561	ENDIF
562	end do
563	c
564	ierr = NF_INQ_VARID (nid, "albe", nvarid)
565	IF (ierr .NE. NF_NOERR) THEN
566	PRINT*, "lect_start_archive: Le champ <albe> est absent"
567	CALL abort
568	ENDIF
569	#ifdef NC_DOUBLE
570	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,albeold)
571	#else
572	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,albeold)
573	#endif
574	IF (ierr .NE. NF_NOERR) THEN
575	PRINT*, "lect_start_archive: Lecture echouee pour <albe>"
576	CALL abort
577	ENDIF
578	c
579	ierr = NF_INQ_VARID (nid, "radsol", nvarid)
580	IF (ierr .NE. NF_NOERR) THEN
581	PRINT*, "lect_start_archive: Le champ <radsol> est absent"
582	CALL abort
583	ENDIF
584	#ifdef NC_DOUBLE
585	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,radsolold)
586	#else
587	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,radsolold)
588	#endif
589	IF (ierr .NE. NF_NOERR) THEN
590	PRINT*, "lect_start_archive: Lecture echouee pour <radsol>"
591	CALL abort
592	ENDIF
593	c
594	ierr = NF_INQ_VARID (nid, "sollw", nvarid)
595	IF (ierr .NE. NF_NOERR) THEN
596	PRINT*, "lect_start_archive: Le champ <sollw> est absent"
597	CALL abort
598	ENDIF
599	#ifdef NC_DOUBLE
600	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,sollwold)
601	#else
602	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,sollwold)
603	#endif
604	IF (ierr .NE. NF_NOERR) THEN
605	PRINT*, "lect_start_archive: Lecture echouee pour <sollw>"
606	CALL abort
607	ENDIF
608	c
609	ierr = NF_INQ_VARID (nid, "solsw", nvarid)
610	IF (ierr .NE. NF_NOERR) THEN
611	PRINT*, "lect_start_archive: Le champ <solsw> est absent"
612	CALL abort
613	ENDIF
614	#ifdef NC_DOUBLE
615	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,solswold)
616	#else
617	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,solswold)
618	#endif
619	IF (ierr .NE. NF_NOERR) THEN
620	PRINT*, "lect_start_archive: Lecture echouee pour <solsw>"
621	CALL abort
622	ENDIF
623	c
624	ierr = NF_INQ_VARID (nid, "fder", nvarid)
625	IF (ierr .NE. NF_NOERR) THEN
626	PRINT*, "lect_start_archive: Le champ <fder> est absent"
627	CALL abort
628	ENDIF
629	#ifdef NC_DOUBLE
630	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,fderold)
631	#else
632	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,fderold)
633	#endif
634	IF (ierr .NE. NF_NOERR) THEN
635	PRINT*, "lect_start_archive: Lecture echouee pour <fder>"
636	CALL abort
637	ENDIF
638	c
639	ierr = NF_INQ_VARID (nid, "dlw", nvarid)
640	IF (ierr .NE. NF_NOERR) THEN
641	PRINT*, "lect_start_archive: Le champ <dlw> est absent"
642	CALL abort
643	ENDIF
644	#ifdef NC_DOUBLE
645	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,dlwold)
646	#else
647	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,dlwold)
648	#endif
649	IF (ierr .NE. NF_NOERR) THEN
650	PRINT*, "lect_start_archive: Lecture echouee pour <dlw>"
651	CALL abort
652	ENDIF
653	c
654	ierr = NF_INQ_VARID (nid, "sollwdown", nvarid)
655	IF (ierr .NE. NF_NOERR) THEN
656	PRINT*, "lect_start_archive: Le champ <sollwdown> est absent"
657	CALL abort
658	ENDIF
659	#ifdef NC_DOUBLE
660	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,sollwdownold)
661	#else
662	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,sollwdownold)
663	#endif
664	IF (ierr .NE. NF_NOERR) THEN
665	PRINT*, "lect_start_archive: Lecture echouee pour <sollwdown>"
666	CALL abort
667	ENDIF
668	c
669	ierr = NF_INQ_VARID (nid, "zmea", nvarid)
670	IF (ierr .NE. NF_NOERR) THEN
671	PRINT*, "lect_start_archive: Le champ <zmea> est absent"
672	PRINT*, " Il est donc initialise a zero"
673	zmeaold=0.
674	ENDIF
675	#ifdef NC_DOUBLE
676	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zmeaold)
677	#else
678	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zmeaold)
679	#endif
680	IF (ierr .NE. NF_NOERR) THEN
681	PRINT*, "lect_start_archive: Lecture echouee pour <zmea>"
682	CALL abort
683	ENDIF
684	c
685	ierr = NF_INQ_VARID (nid, "zstd", nvarid)
686	IF (ierr .NE. NF_NOERR) THEN
687	PRINT*, "lect_start_archive: Le champ <zstd> est absent"
688	PRINT*, " Il est donc initialise a zero"
689	zstdold=0.
690	ENDIF
691	#ifdef NC_DOUBLE
692	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zstdold)
693	#else
694	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zstdold)
695	#endif
696	IF (ierr .NE. NF_NOERR) THEN
697	PRINT*, "lect_start_archive: Lecture echouee pour <zstd>"
698	CALL abort
699	ENDIF
700	c
701	ierr = NF_INQ_VARID (nid, "zsig", nvarid)
702	IF (ierr .NE. NF_NOERR) THEN
703	PRINT*, "lect_start_archive: Le champ <zsig> est absent"
704	PRINT*, " Il est donc initialise a zero"
705	zsigold=0.
706	ENDIF
707	#ifdef NC_DOUBLE
708	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zsigold)
709	#else
710	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zsigold)
711	#endif
712	IF (ierr .NE. NF_NOERR) THEN
713	PRINT*, "lect_start_archive: Lecture echouee pour <zsig>"
714	CALL abort
715	ENDIF
716	c
717	ierr = NF_INQ_VARID (nid, "zgam", nvarid)
718	IF (ierr .NE. NF_NOERR) THEN
719	PRINT*, "lect_start_archive: Le champ <zgam> est absent"
720	PRINT*, " Il est donc initialise a zero"
721	zgamold=0.
722	ENDIF
723	#ifdef NC_DOUBLE
724	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zgamold)
725	#else
726	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zgamold)
727	#endif
728	IF (ierr .NE. NF_NOERR) THEN
729	PRINT*, "lect_start_archive: Lecture echouee pour <zgam>"
730	CALL abort
731	ENDIF
732	c
733	ierr = NF_INQ_VARID (nid, "zthe", nvarid)
734	IF (ierr .NE. NF_NOERR) THEN
735	PRINT*, "lect_start_archive: Le champ <zthe> est absent"
736	PRINT*, " Il est donc initialise a zero"
737	ztheold=0.
738	ENDIF
739	#ifdef NC_DOUBLE
740	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,ztheold)
741	#else
742	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,ztheold)
743	#endif
744	IF (ierr .NE. NF_NOERR) THEN
745	PRINT*, "lect_start_archive: Lecture echouee pour <zthe>"
746	CALL abort
747	ENDIF
748	c
749	ierr = NF_INQ_VARID (nid, "zpic", nvarid)
750	IF (ierr .NE. NF_NOERR) THEN
751	PRINT*, "lect_start_archive: Le champ <zpic> est absent"
752	PRINT*, " Il est donc initialise a zero"
753	zpicold=0.
754	ENDIF
755	#ifdef NC_DOUBLE
756	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zpicold)
757	#else
758	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zpicold)
759	#endif
760	IF (ierr .NE. NF_NOERR) THEN
761	PRINT*, "lect_start_archive: Lecture echouee pour <zpic>"
762	CALL abort
763	ENDIF
764	c
765	ierr = NF_INQ_VARID (nid, "zval", nvarid)
766	IF (ierr .NE. NF_NOERR) THEN
767	PRINT*, "lect_start_archive: Le champ <zval> est absent"
768	PRINT*, " Il est donc initialise a zero"
769	zvalold=0.
770	ENDIF
771	#ifdef NC_DOUBLE
772	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,zvalold)
773	#else
774	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,zvalold)
775	#endif
776	IF (ierr .NE. NF_NOERR) THEN
777	PRINT*, "lect_start_archive: Lecture echouee pour <zval>"
778	CALL abort
779	ENDIF
780	c
781
782	print*,"Lecture champs 2D OK"
783
784	c-----------------------------------------------------------------------
785	c Lecture des champs 3D ()
786	c-----------------------------------------------------------------------
787
788	start=(/1,1,1,1/)
789	counter=(/imold+1,jmold+1,lmold,1/)
790	c
791	ierr = NF_INQ_VARID (nid,"temp", nvarid)
792	IF (ierr .NE. NF_NOERR) THEN
793	PRINT*, "lect_start_archive: Le champ <temp> est absent"
794	CALL abort
795	ENDIF
796	#ifdef NC_DOUBLE
797	ierr = NF_GET_VARA_DOUBLE(nid, nvarid, start, counter, told)
798	#else
799	ierr = NF_GET_VARA_REAL(nid, nvarid, start, counter, told)
800	#endif
801	IF (ierr .NE. NF_NOERR) THEN
802	PRINT*, "lect_start_archive: Lecture echouee pour <temp>"
803	CALL abort
804	ENDIF
805	c
806	ierr = NF_INQ_VARID (nid,"u", nvarid)
807	IF (ierr .NE. NF_NOERR) THEN
808	PRINT*, "lect_start_archive: Le champ <u> est absent"
809	CALL abort
810	ENDIF
811	#ifdef NC_DOUBLE
812	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,uold)
813	#else
814	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,uold)
815	#endif
816	IF (ierr .NE. NF_NOERR) THEN
817	PRINT*, "lect_start_archive: Lecture echouee pour <u>"
818	CALL abort
819	ENDIF
820	c
821	ierr = NF_INQ_VARID (nid,"v", nvarid)
822	IF (ierr .NE. NF_NOERR) THEN
823	PRINT*, "lect_start_archive: Le champ <v> est absent"
824	CALL abort
825	ENDIF
826	#ifdef NC_DOUBLE
827	ierr = NF_GET_VARA_DOUBLE(nid, nvarid,start,counter,vold)
828	#else
829	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,vold)
830	#endif
831	IF (ierr .NE. NF_NOERR) THEN
832	PRINT*, "lect_start_archive: Lecture echouee pour <v>"
833	CALL abort
834	ENDIF
835	c
836
837	c TNAME: IL EST LU A PARTIR DE traceur.def (mettre l'ancien si
838	c changement du nombre de traceurs)
839
840	IF(iflag_trac.eq.1) THEN
841	DO iq=1,nqtot
842	ierr = NF_INQ_VARID (nid, tname(iq), nvarid)
843	IF (ierr .NE. NF_NOERR) THEN
844	PRINT*, "dynetat0: Le champ <"//tname(iq)//"> est absent"
845	PRINT*, " Il est donc initialise a zero"
846	qold(:,:,:,iq)=0.
847	ELSE
848	#ifdef NC_DOUBLE
849	ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,counter,qold(1,1,1,iq))
850	#else
851	ierr = NF_GET_VARA_REAL(nid, nvarid,start,counter,qold(1,1,1,iq))
852	#endif
853	IF (ierr .NE. NF_NOERR) THEN
854	PRINT*, "dynetat0: Lecture echouee pour "//tname(iq)
855	CALL abort
856	ENDIF
857	ENDIF
858	ENDDO
859	ENDIF
860
861
862	print*,"Lecture champs 3D OK"
863
864	c=======================================================================
865	c INTERPOLATION DANS LA NOUVELLE GRILLE et initialisation des variables
866	c=======================================================================
867	c Interpolation horizontale puis passage dans la grille physique pour
868	c les variables physique
869	c Interpolation verticale puis horizontale pour chaque variable 3D
870	c=======================================================================
871
872	c-----------------------------------------------------------------------
873	c Variable 2d :
874	c-----------------------------------------------------------------------
875
876	c Relief
877	! topoflag = F: we keep the existing topography
878	! = T: we read it from the Relief.nc file
879	! topoflag need to be in the specific run.def for newstart
880	topoflag = . FALSE .
881	CALL getin('topoflag',topoflag)
882
883	print*,zmeaold(2,1:10)
884
885	IF ( topoflag ) then
886	print*,"Topography (phis) read in file Relief.nc"
887	print*,"For Venus, map directly inverted in Relief.nc"
888	CALL startget_phys2d('surfgeo',iip1,jjp1,rlonv,rlatu,phis,0.0,
889	. jjm ,rlonu,rlatv,.true.)
890	c needed ?
891	phis(iip1,:) = phis(1,:)
892
893	CALL startget_phys1d('zmea',iip1,jjp1,rlonv,rlatu,ngridmx,zmea,
894	. 0.0,jjm,rlonu,rlatv,.true.)
895	CALL startget_phys1d('zstd',iip1,jjp1,rlonv,rlatu,ngridmx,zstd,
896	. 0.0,jjm,rlonu,rlatv,.true.)
897	CALL startget_phys1d('zsig',iip1,jjp1,rlonv,rlatu,ngridmx,zsig,
898	. 0.0,jjm,rlonu,rlatv,.true.)
899	CALL startget_phys1d('zgam',iip1,jjp1,rlonv,rlatu,ngridmx,zgam,
900	. 0.0,jjm,rlonu,rlatv,.true.)
901	CALL startget_phys1d('zthe',iip1,jjp1,rlonv,rlatu,ngridmx,zthe,
902	. 0.0,jjm,rlonu,rlatv,.true.)
903	CALL startget_phys1d('zpic',iip1,jjp1,rlonv,rlatu,ngridmx,zpic,
904	. 0.0,jjm,rlonu,rlatv,.true.)
905	CALL startget_phys1d('zval',iip1,jjp1,rlonv,rlatu,ngridmx,zval,
906	. 0.0,jjm,rlonu,rlatv,.true.)
907
908	ELSE
909	print*,'Using existing topography'
910	call interp_horiz (phisold,phis,imold,jmold,iim,jjm,1,
911	& rlonuold,rlatvold,rlonu,rlatv)
912
913	call interp_horiz (zmeaold,zmeaS,imold,jmold,iim,jjm,1,
914	& rlonuold,rlatvold,rlonu,rlatv)
915	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zmeaS,zmea)
916	call interp_horiz (zstdold,zstdS,imold,jmold,iim,jjm,1,
917	& rlonuold,rlatvold,rlonu,rlatv)
918	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zstdS,zstd)
919	call interp_horiz (zsigold,zsigS,imold,jmold,iim,jjm,1,
920	& rlonuold,rlatvold,rlonu,rlatv)
921	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zsigS,zsig)
922	call interp_horiz (zgamold,zgamS,imold,jmold,iim,jjm,1,
923	& rlonuold,rlatvold,rlonu,rlatv)
924	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zgamS,zgam)
925	call interp_horiz (ztheold,ztheS,imold,jmold,iim,jjm,1,
926	& rlonuold,rlatvold,rlonu,rlatv)
927	call gr_dyn_fi (1,iip1,jjp1,ngridmx,ztheS,zthe)
928	call interp_horiz (zpicold,zpicS,imold,jmold,iim,jjm,1,
929	& rlonuold,rlatvold,rlonu,rlatv)
930	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zpicS,zpic)
931	call interp_horiz (zvalold,zvalS,imold,jmold,iim,jjm,1,
932	& rlonuold,rlatvold,rlonu,rlatv)
933	call gr_dyn_fi (1,iip1,jjp1,ngridmx,zvalS,zval)
934	ENDIF
935
936	print*,"New surface geopotential OK"
937
938	c Lat et lon pour physique
939	do i=1,iip1
940	rlatS(i,:)=rlatu(:)*180./pi
941	enddo
942	call gr_dyn_fi (1,iip1,jjp1,ngridmx,rlatS,rlat)
943
944	do j=2,jjm
945	rlonS(:,j)=rlonv(:)*180./pi
946	enddo
947	rlonS(:,1)=0.
948	rlonS(:,jjp1)=0.
949	call gr_dyn_fi (1,iip1,jjp1,ngridmx,rlonS,rlon)
950
951	c Temperature de surface
952	call interp_horiz (tsurfold,tsurfS,imold,jmold,iim,jjm,1,
953	& rlonuold,rlatvold,rlonu,rlatv)
954	call gr_dyn_fi (1,iip1,jjp1,ngridmx,tsurfS,tsurf)
955	c write(44,*) 'tsurf', tsurf
956
957	c Temperature du sous-sol
958	call interp_horiz(tsoilold,tsoilS,
959	& imold,jmold,iim,jjm,nsoilmx,
960	& rlonuold,rlatvold,rlonu,rlatv)
961	call gr_dyn_fi (nsoilmx,iip1,jjp1,ngridmx,tsoilS,tsoil)
962	c write(45,*) 'tsoil',tsoil
963
964	! CHANGING ALBEDO: may be done through run.def
965	albedoflag = . FALSE .
966	CALL getin('albedoflag',albedoflag)
967
968	IF ( albedoflag ) then
969	print*,"Albedo is reinitialized to the albedo value in run.def"
970	print*,"We may want to consider a map later on..."
971	albedo=0.1
972	CALL getin('albedo',albedo)
973	albe=albedo
974	ELSE
975	call interp_horiz (albeold,albeS,imold,jmold,iim,jjm,1,
976	& rlonuold,rlatvold,rlonu,rlatv)
977	call gr_dyn_fi (1,iip1,jjp1,ngridmx,albeS,albe)
978	ENDIF
979
980	call interp_horiz (radsolold,radsolS,imold,jmold,iim,jjm,1,
981	& rlonuold,rlatvold,rlonu,rlatv)
982	call gr_dyn_fi (1,iip1,jjp1,ngridmx,radsolS,radsol)
983
984	call interp_horiz (sollwold,sollwS,imold,jmold,iim,jjm,1,
985	& rlonuold,rlatvold,rlonu,rlatv)
986	call gr_dyn_fi (1,iip1,jjp1,ngridmx,sollwS,sollw)
987
988	call interp_horiz (solswold,solswS,imold,jmold,iim,jjm,1,
989	& rlonuold,rlatvold,rlonu,rlatv)
990	call gr_dyn_fi (1,iip1,jjp1,ngridmx,solswS,solsw)
991
992	call interp_horiz (fderold,fderS,imold,jmold,iim,jjm,1,
993	& rlonuold,rlatvold,rlonu,rlatv)
994	call gr_dyn_fi (1,iip1,jjp1,ngridmx,fderS,fder)
995
996	call interp_horiz (dlwold,dlwS,imold,jmold,iim,jjm,1,
997	& rlonuold,rlatvold,rlonu,rlatv)
998	call gr_dyn_fi (1,iip1,jjp1,ngridmx,dlwS,dlw)
999
1000	call interp_horiz (sollwdownold,sollwdownS,imold,jmold,iim,jjm,1,
1001	& rlonuold,rlatvold,rlonu,rlatv)
1002	call gr_dyn_fi (1,iip1,jjp1,ngridmx,sollwdownS,sollwdown)
1003
1004	print*,"Nouvelles var physiques OK"
1005
1006	c-----------------------------------------------------------------------
1007	c Traitement special de la pression au sol :
1008	c-----------------------------------------------------------------------
1009
1010	c Extrapolation la pression dans la nouvelle grille
1011	call interp_horiz(psold,ps,imold,jmold,iim,jjm,1,
1012	& rlonuold,rlatvold,rlonu,rlatv)
1013
1014	c On assure la conservation de la masse de l'atmosphere
1015	c--------------------------------------------------------------
1016
1017	ptotal = 0.
1018	DO j=1,jjp1
1019	DO i=1,iim
1020	ptotal=ptotal+ps(i,j)*aire(i,j)/g
1021	END DO
1022	END DO
1023
1024	write(,)
1025	write(,)'Ancienne grille: masse de l atm :',ptotalold
1026	write(,)'Nouvelle grille: masse de l atm :',ptotal
1027	write (,) 'Ratio new atm./ old atm =', ptotal/ptotalold
1028	write(,)
1029
1030
1031	DO j=1,jjp1
1032	DO i=1,iip1
1033	ps(i,j)=ps(i,j) * ptotalold/ptotal
1034	END DO
1035	END DO
1036
1037	c la pression de surface et les temperatures ne sont pas reequilibrees en fonction
1038	c de la nouvelle topographie...
1039	c Si l'ajustement inevitable du debut pose des problemes, voir le newstart martien.
1040
1041	print*,"Nouvelle ps OK"
1042	print*,"If changes done on topography, beware !"
1043	print*,"Some time may be needed for adjustments at the beginning"
1044	print*,"so if unstable, relax the timestep and/or dissipation."
1045
1046	c-----------------------------------------------------------------------
1047	c Variable 3d :
1048	c-----------------------------------------------------------------------
1049
1050	CALL pression(ip1jmp1, ap, bp, ps, p3d)
1051	if (disvert_type==1) then
1052	CALL exner_hyb( ip1jmp1, ps, p3d, pks, pk, pkf )
1053	else ! we assume that we are in the disvert_type==2 case
1054	CALL exner_milieu( ip1jmp1, ps, p3d, pks, pk, pkf )
1055	endif
1056
1057	c temperatures atmospheriques
1058
1059	c ATTENTION: peut servir, mais bon...
1060	c modif: profil uniforme
1061	c do l=1,lmold
1062	c do j=1,jmold+1
1063	c do i=1,imold+1
1064	c told(i,j,l)=told(1,jmold/2,l)
1065	c enddo
1066	c enddo
1067	c enddo
1068
1069	write (,) 'told ', told (1,jmold+1,1) ! INFO
1070	call interp_vert
1071	& (told,var,lmold,llm,apold,bpold,ap,bp,
1072	& psold,(imold+1)*(jmold+1))
1073	write (,) 'var ', var (1,jmold+1,1) ! INFO
1074	call interp_horiz(var,t,imold,jmold,iim,jjm,llm,
1075	& rlonuold,rlatvold,rlonu,rlatv)
1076	write (,) 'T ', T(1,jjp1,1) ! INFO
1077	! pour info:
1078	! Si extension verticale, la T est extrapolee constante au-dessus de lmold
1079
1080	c passage grille physique pour restartphy.nc
1081	call gr_dyn_fi (llm,iip1,jjp1,ngridmx,T,t_fi)
1082
1083	! ADAPTATION GCM POUR CP(T)
1084	c passage en temperature potentielle
1085	call t2tpot(ip1jmp1*llm,T,teta,pk)
1086	c on assure la periodicite
1087	teta(iip1,:,:) = teta(1,:,:)
1088
1089	! RESETING U TO 0: may be done through run.def
1090	razvitu = . FALSE .
1091	CALL getin('razvitu',razvitu)
1092	razvitv = . FALSE .
1093	CALL getin('razvitv',razvitv)
1094
1095	c calcul des champ de vent; passage en vent covariant
1096	write (,) 'uold ', uold (1,2,1) ! INFO
1097	call interp_vert
1098	& (uold,var,lmold,llm,apold,bpold,ap,bp,
1099	& psold,(imold+1)*(jmold+1))
1100	write (,) 'var ', var (1,2,1) ! INFO
1101	call interp_horiz(var,us,imold,jmold,iim,jjm,llm,
1102	& rlonuold,rlatvold,rlonu,rlatv)
1103	write (,) 'us ', us (1,2,1) ! INFO
1104
1105	call interp_vert
1106	& (vold,var,lmold,llm,
1107	& apold,bpold,ap,bp,psold,(imold+1)*(jmold+1))
1108	call interp_horiz(var,vs,imold,jmold,iim,jjm,llm,
1109	& rlonuold,rlatvold,rlonu,rlatv)
1110	call scal_wind(us,vs,unat,vnat)
1111	! Reseting u=0
1112	if (razvitu) then
1113	unat(:,:,:) = 0.
1114	endif
1115	write (,) 'unat ', unat (1,2,1) ! INFO
1116	do l=1,llm
1117	do j = 1, jjp1
1118	do i=1,iip1
1119	ucov( i,j,l ) = unat( i,j,l ) * cu(i,j)
1120	! pour info:
1121	! Si extension verticale, on impose u=0 au-dessus de lmold
1122	if (l.gt.lmold) ucov( i,j,l ) = 0
1123	end do
1124	end do
1125	end do
1126	write (,) 'ucov ', ucov (1,2,1) ! INFO
1127	c write(48,*) 'ucov',ucov
1128	! Reseting v=0
1129	if (razvitv) then
1130	vnat(:,:,:) = 0.
1131	endif
1132	write (,) 'vnat ', vnat (1,2,1) ! INFO
1133	do l=1,llm
1134	do j = 1, jjm
1135	do i=1,iim
1136	vcov( i,j,l ) = vnat( i,j,l ) * cv(i,j)
1137	! pour info:
1138	! Si extension verticale, on impose v=0 au-dessus de lmold
1139	if (l.gt.lmold) vcov( i,j,l ) = 0
1140	end do
1141	vcov( iip1,j,l ) = vcov( 1,j,l )
1142	end do
1143	end do
1144	c write(49,*) 'ucov',vcov
1145
1146	c masse: on la recalcule (ps a été ajustée pour conserver la masse totale)
1147	call massdair(p3d,masse)
1148
1149	c traceurs 3D
1150	do iq = 1, nqtot
1151	call interp_vert(qold(1,1,1,iq),var,lmold,llm,
1152	& apold,bpold,ap,bp,psold,(imold+1)*(jmold+1))
1153	call interp_horiz(var,q(1,1,1,iq),imold,jmold,iim,jjm,llm,
1154	& rlonuold,rlatvold,rlonu,rlatv)
1155	enddo
1156
1157	print*,"Nouvelles var dynamiques OK"
1158
1159	c=======================================================================
1160	c Ecriture des restart.nc et restartphy.nc :
1161	c=======================================================================
1162
1163	call writerestart('restart.nc',tab_cntrl_dyn,
1164	. phis,vcov,ucov,teta,masse,q,ps)
1165
1166	print*,"restart OK"
1167
1168	call writerestartphy('restartphy.nc',tab_cntrl_fi,ngridmx,llm,
1169	. rlat,rlon,tsurf,tsoil,albe,solsw, sollw,fder,dlw,
1170	. sollwdown,radsol,
1171	. zmea,zstd,zsig,zgam,zthe,zpic,zval,
1172	. t_fi)
1173
1174	print*,"restartphy OK"
1175
1176	end

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format