Context Navigation

newstart.F @ 937

Last change on this file since 937 was 927, checked in by slebonnois, 12 years ago
SL: mise a jour de quelques docs + correction bug dans newstart (Venus)
File size: 37.1 KB

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

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

Download in other formats:

Original Format