Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

startvar.F @ 5416

Last change on this file since 5416 was 253, checked in by (none), 23 years ago
This commit was manufactured by cvs2svn to create branch 'rel-1-0-patch'.
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 32.6 KB

Rev	Line
[2]	1	MODULE startvar
	2	!
	3	!
	4	! There are three ways to access data from the database of atmospheric data which
	5	! can be used to initialize the model. This depends on the type of field which needs
	6	! to be extracted. In any case the call should come after a restget and should be of the type :
	7	! CALL startget(...)
	8	!
	9	! We will details the possible arguments to startget here :
	10	!
	11	! - A 2D variable on the dynamical grid :
	12	! CALL startget(varname, iml, jml, lon_in, lat_in, champ, val_ex)
	13	!
	14	!
	15	! - A 1D variable on the physical grid :
	16	! CALL startget(varname, iml, jml, lon_in, lat_in, nbindex, champ, val_exp)
	17	!
	18	!
	19	! - A 3D variable on the dynamical grid :
	20	! CALL startget(varname, iml, jml, lon_in, lat_in, lml, pls, workvar, champ, val_exp)
	21	!
	22	!
	23	! There is special constraint on the atmospheric data base except that the
	24	! the data needs to be in netCDF and the variables should have the the following
	25	! names in the file :
	26	!
	27	! 'RELIEF' : High resolution orography
	28	! 'ST' : Surface temperature
	29	! 'CDSW' : Soil moisture
	30	! 'Z' : Surface geopotential
	31	! 'SP' : Surface pressure
	32	! 'U' : East ward wind
	33	! 'V' : Northward wind
	34	! 'TEMP' : Temperature
	35	! 'R' : Relative humidity
	36	!
	37	USE ioipsl
	38	!
	39	!
	40	IMPLICIT NONE
	41	!
	42	!
	43	PRIVATE
	44	PUBLIC startget
	45	!
	46	!
	47	INTERFACE startget
	48	MODULE PROCEDURE startget_phys2d, startget_phys1d, startget_dyn
	49	END INTERFACE
	50	!
	51	INTEGER, SAVE :: fid_phys, fid_dyn
	52	INTEGER, SAVE :: iml_phys, iml_rel, iml_dyn
	53	INTEGER, SAVE :: jml_phys, jml_rel, jml_dyn
	54	INTEGER, SAVE :: llm_dyn, ttm_dyn
	55	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: lon_phys, lon_rug,
	56	. lon_alb, lon_rel, lon_dyn
	57	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: lat_phys, lat_rug,
	58	. lat_alb, lat_rel, lat_dyn
	59	REAL, ALLOCATABLE, SAVE, DIMENSION (:) :: lev_dyn
	60	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: relief, zstd, zsig,
	61	. zgam, zthe, zpic, zval
	62	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: rugo, masque, phis
	63	!
	64	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:) :: tsol, qsol, psol_dyn
	65	REAL, ALLOCATABLE, SAVE, DIMENSION (:,:,:) :: var_ana3d
	66	!
	67	CONTAINS
	68	!
	69	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	70	!
	71	SUBROUTINE startget_phys2d(varname, iml, jml, lon_in, lat_in,
	72	. champ, val_exp)
	73	!
	74	! There is a big mess with the size in logitude, should it be iml or iml+1.
	75	! I have chosen to use the iml+1 as an argument to this routine and we declare
	76	! internaly smaler fields when needed. This needs to be cleared once and for all in LMDZ.
	77	! A convention is required.
	78	!
	79	!
	80	CHARACTER(), INTENT(in) :: varname
	81	INTEGER, INTENT(in) :: iml, jml
	82	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	83	REAL, INTENT(inout) :: champ(iml,jml)
	84	REAL, INTENT(in) :: val_exp
	85	!
	86	! This routine only works if the variable does not exist or is constant
	87	!
	88	IF ( MINVAL(champ(:,:)).EQ.MAXVAL(champ(:,:)) .AND.
	89	.MINVAL(champ(:,:)).EQ.val_exp ) THEN
	90	!
	91	SELECTCASE(varname)
	92	!
	93	CASE ('relief')
	94	!
	95	! If we do not have the orography we need to get it
	96	!
	97	IF ( .NOT.ALLOCATED(relief)) THEN
	98	!
	99	CALL start_init_orog( iml, jml, lon_in, lat_in)
	100	!
	101	ENDIF
	102	!
	103	IF (SIZE(relief) .NE. SIZE(champ)) THEN
	104	!
	105	WRITE(,) 'STARTVAR module has been',
	106	.' initialized to the wrong size'
	107	STOP
	108	!
	109	ENDIF
	110	!
	111	champ(:,:) = relief(:,:)
	112	!
	113	CASE ('rugosite')
	114	!
	115	! If we do not have the orography we need to get it
	116	!
	117	IF ( .NOT.ALLOCATED(rugo)) THEN
	118	!
	119	CALL start_init_orog( iml, jml, lon_in, lat_in)
	120	!
	121	ENDIF
	122	!
	123	IF (SIZE(rugo) .NE. SIZE(champ)) THEN
	124	!
	125	WRITE(,)
	126	. 'STARTVAR module has been initialized to the wrong size'
	127	STOP
	128	!
	129	ENDIF
	130	!
	131	champ(:,:) = rugo(:,:)
	132	!
	133	CASE ('masque')
	134	!
	135	! If we do not have the orography we need to get it
	136	!
	137	IF ( .NOT.ALLOCATED(masque)) THEN
	138	!
	139	CALL start_init_orog( iml, jml, lon_in, lat_in)
	140	!
	141	ENDIF
	142	!
	143	IF (SIZE(masque) .NE. SIZE(champ)) THEN
	144	!
	145	WRITE(,)
	146	. 'STARTVAR module has been initialized to the wrong size'
	147	STOP
	148	!
	149	ENDIF
	150	!
	151	champ(:,:) = masque(:,:)
	152	!
	153	CASE ('surfgeo')
	154	!
	155	! If we do not have the orography we need to get it
	156	!
	157	IF ( .NOT.ALLOCATED(phis)) THEN
	158	!
	159	CALL start_init_orog( iml, jml, lon_in, lat_in)
	160	!
	161	ENDIF
	162	!
	163	IF (SIZE(phis) .NE. SIZE(champ)) THEN
	164	!
	165	WRITE(,)
	166	. 'STARTVAR module has been initialized to the wrong size'
	167	STOP
	168	!
	169	ENDIF
	170	!
	171	champ(:,:) = phis(:,:)
	172	!
	173	CASE ('psol')
	174	!
	175	! If we do not have the orography we need to get it
	176	!
	177	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	178	!
	179	CALL start_init_dyn( iml, jml, lon_in, lat_in)
	180	!
	181	ENDIF
	182	!
	183	IF (SIZE(psol_dyn) .NE. SIZE(champ)) THEN
	184	!
	185	WRITE(,)
	186	. 'STARTVAR module has been initialized to the wrong size'
	187	STOP
	188	!
	189	ENDIF
	190	!
	191	champ(:,:) = psol_dyn(:,:)
	192	!
	193	CASE DEFAULT
	194	!
	195	WRITE(,) 'startget_phys2d'
	196	WRITE(,) 'No rule is present to extract variable',
	197	. varname(:LEN_TRIM(varname)),' from any data set'
	198	STOP
	199	!
	200	END SELECT
	201	!
	202	ELSE
	203	!
	204	! There are a few fields we might need if we need to interpolate 3D filed. Thus if they come through here we
	205	! will catch them
	206	!
	207	SELECTCASE(varname)
	208	!
	209	CASE ('surfgeo')
	210	!
	211	IF ( .NOT.ALLOCATED(phis)) THEN
	212	ALLOCATE(phis(iml,jml))
	213	ENDIF
	214	!
	215	IF (SIZE(phis) .NE. SIZE(champ)) THEN
	216	!
	217	WRITE(,)
	218	. 'STARTVAR module has been initialized to the wrong size'
	219	STOP
	220	!
	221	ENDIF
	222	!
	223	phis(:,:) = champ(:,:)
	224	!
	225	END SELECT
	226	!
	227	ENDIF
	228	!
	229	END SUBROUTINE startget_phys2d
	230	!
	231	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	232	!
	233	SUBROUTINE start_init_orog( iml, jml, lon_in, lat_in)
	234	!
	235	INTEGER, INTENT(in) :: iml, jml
	236	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	237	!
	238	! LOCAL
	239	!
	240	REAL :: lev(1), date, dt
	241	INTEGER :: itau(1), fid
	242	INTEGER :: llm_tmp, ttm_tmp
	243	INTEGER :: i, j
	244	INTEGER :: iret
	245	REAL, ALLOCATABLE :: relief_hi(:,:)
	246	REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
	247	REAL, ALLOCATABLE :: tmp_var(:,:)
	248	INTEGER, ALLOCATABLE :: tmp_int(:,:)
	249	!
	250	CHARACTER*120 :: orogfname
	251	LOGICAL :: check=.TRUE.
	252	!
	253	!
	254	orogfname = 'Relief.nc'
	255	!
	256	IF ( check ) WRITE(,) 'Reading the high resolution orography'
	257	!
	258	CALL flininfo(orogfname,iml_rel, jml_rel, llm_tmp, ttm_tmp, fid)
	259	!
	260	ALLOCATE (lat_rel(iml_rel,jml_rel), stat=iret)
	261	ALLOCATE (lon_rel(iml_rel,jml_rel), stat=iret)
	262	ALLOCATE (relief_hi(iml_rel,jml_rel), stat=iret)
	263	!
	264	CALL flinopen(orogfname, .FALSE., iml_rel, jml_rel,
	265	.llm_tmp, lon_rel, lat_rel, lev, ttm_tmp,
	266	. itau, date, dt, fid)
	267	!
	268	CALL flinget(fid, 'RELIEF', iml_rel, jml_rel, llm_tmp,
	269	. ttm_tmp, 1, 1, relief_hi)
	270	!
	271	CALL flinclo(fid)
	272	!
	273	! In case we have a file which is in degrees we do the transformation
	274	!
	275	ALLOCATE(lon_rad(iml_rel))
	276	IF ( MAXVAL(lon_rel(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	277	lon_rad(:) = lon_rel(:,1) * 2.0 * ASIN(1.0) / 180.0
	278	ELSE
	279	lon_rad(:) = lon_rel(:,1)
	280	ENDIF
	281	ALLOCATE(lat_rad(jml_rel))
	282	IF ( MAXVAL(lat_rel(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	283	lat_rad(:) = lat_rel(1,:) * 2.0 * ASIN(1.0) / 180.0
	284	ELSE
	285	lat_rad(:) = lat_rel(1,:)
	286	ENDIF
	287	!
	288	!
	289	IF ( check ) WRITE(,) 'Computes all the parameters needed',
	290	.' for the gravity wave drag code'
	291	!
	292	! Allocate the data we need to put in the interpolated fields
	293	!
	294	! RELIEF: orographie moyenne
	295	ALLOCATE(relief(iml,jml))
	296	! zphi : orographie moyenne
	297	ALLOCATE(phis(iml,jml))
	298	! zstd: deviation standard de l'orographie sous-maille
	299	ALLOCATE(zstd(iml,jml))
	300	! zsig: pente de l'orographie sous-maille
	301	ALLOCATE(zsig(iml,jml))
	302	! zgam: anisotropy de l'orographie sous maille
	303	ALLOCATE(zgam(iml,jml))
	304	! zthe: orientation de l'axe oriente dans la direction
	305	! de plus grande pente de l'orographie sous maille
	306	ALLOCATE(zthe(iml,jml))
	307	! zpic: hauteur pics de la SSO
	308	ALLOCATE(zpic(iml,jml))
	309	! zval: hauteur vallees de la SSO
	310	ALLOCATE(zval(iml,jml))
	311	! masque : Masque terre ocean
	312	ALLOCATE(tmp_int(iml,jml))
	313	ALLOCATE(masque(iml,jml))
	314	!
	315	CALL grid_noro(iml_rel, jml_rel, lon_rad, lat_rad, relief_hi,
	316	. iml-1, jml, lon_in, lat_in,
	317	. phis, relief, zstd, zsig, zgam, zthe, zpic, zval, tmp_int)
[24]	318	phis = phis * 9.81
[2]	319	!
	320	masque(:,:) = FLOAT(tmp_int(:,:))
	321	!
	322	! Compute surface roughness
	323	!
	324	IF ( check ) WRITE(,)
	325	.'Compute surface roughness induced by the orography'
	326	!
	327	ALLOCATE(rugo(iml,jml))
	328	ALLOCATE(tmp_var(iml-1,jml))
	329	!
	330	CALL rugsoro(iml_rel, jml_rel, lon_rad, lat_rad, relief_hi,
	331	. iml-1, jml, lon_in, lat_in, tmp_var)
	332	!
	333	DO j = 1, jml
	334	DO i = 1, iml-1
	335	rugo(i,j) = tmp_var(i,j)
	336	ENDDO
	337	rugo(iml,j) = tmp_var(1,j)
	338	ENDDO
	339	!
	340	! Build land-sea mask
	341	!
	342	!
	343	RETURN
	344	!
	345	END SUBROUTINE start_init_orog
	346	!
	347	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	348	!
	349	SUBROUTINE startget_phys1d(varname, iml, jml, lon_in,
	350	.lat_in, nbindex, champ, val_exp)
	351	!
	352	CHARACTER(), INTENT(in) :: varname
	353	INTEGER, INTENT(in) :: iml, jml, nbindex
	354	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	355	REAL, INTENT(inout) :: champ(nbindex)
	356	REAL, INTENT(in) :: val_exp
	357	!
	358	!
	359	! This routine only works if the variable does not exist or is constant
	360	!
	361	IF ( MINVAL(champ(:)).EQ.MAXVAL(champ(:)) .AND.
	362	.MINVAL(champ(:)).EQ.val_exp ) THEN
	363	SELECTCASE(varname)
	364	CASE ('tsol')
	365	IF ( .NOT.ALLOCATED(tsol)) THEN
	366	CALL start_init_phys( iml, jml, lon_in, lat_in)
	367	ENDIF
	368	IF ( SIZE(tsol) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	369	WRITE(,)
	370	. 'STARTVAR module has been initialized to the wrong size'
	371	STOP
	372	ENDIF
	373	CALL gr_dyn_fi(1, iml, jml, nbindex, tsol, champ)
	374	CASE ('qsol')
	375	IF ( .NOT.ALLOCATED(qsol)) THEN
	376	CALL start_init_phys( iml, jml, lon_in, lat_in)
	377	ENDIF
	378	IF ( SIZE(qsol) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	379	WRITE(,)
	380	. 'STARTVAR module has been initialized to the wrong size'
	381	STOP
	382	ENDIF
	383	CALL gr_dyn_fi(1, iml, jml, nbindex, qsol, champ)
	384	CASE ('psol')
	385	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	386	CALL start_init_dyn( iml, jml, lon_in, lat_in)
	387	ENDIF
	388	IF (SIZE(psol_dyn) .NE. SIZE(lon_in)*SIZE(lat_in)) THEN
	389	WRITE(,)
	390	. 'STARTVAR module has been initialized to the wrong size'
	391	STOP
	392	ENDIF
	393	CALL gr_dyn_fi(1, iml, jml, nbindex, psol_dyn, champ)
	394	CASE ('zmea')
	395	IF ( .NOT.ALLOCATED(relief)) THEN
	396	CALL start_init_orog( iml, jml, lon_in, lat_in)
	397	ENDIF
	398	IF ( SIZE(relief) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	399	WRITE(,)
	400	. 'STARTVAR module has been initialized to the wrong size'
	401	STOP
	402	ENDIF
	403	CALL gr_dyn_fi(1, iml, jml, nbindex, relief, champ)
	404	CASE ('zstd')
	405	IF ( .NOT.ALLOCATED(zstd)) THEN
	406	CALL start_init_orog( iml, jml, lon_in, lat_in)
	407	ENDIF
	408	IF ( SIZE(zstd) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	409	WRITE(,)
	410	. 'STARTVAR module has been initialized to the wrong size'
	411	STOP
	412	ENDIF
	413	CALL gr_dyn_fi(1, iml, jml, nbindex,zstd, champ)
	414	CASE ('zsig')
	415	IF ( .NOT.ALLOCATED(zsig)) THEN
	416	CALL start_init_orog( iml, jml, lon_in, lat_in)
	417	ENDIF
	418	IF ( SIZE(zsig) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	419	WRITE(,)
	420	. 'STARTVAR module has been initialized to the wrong size'
	421	STOP
	422	ENDIF
	423	CALL gr_dyn_fi(1, iml, jml, nbindex,zsig, champ)
	424	CASE ('zgam')
	425	IF ( .NOT.ALLOCATED(zgam)) THEN
	426	CALL start_init_orog( iml, jml, lon_in, lat_in)
	427	ENDIF
	428	IF ( SIZE(zgam) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	429	WRITE(,)
	430	. 'STARTVAR module has been initialized to the wrong size'
	431	STOP
	432	ENDIF
	433	CALL gr_dyn_fi(1, iml, jml, nbindex,zgam, champ)
	434	CASE ('zthe')
	435	IF ( .NOT.ALLOCATED(zthe)) THEN
	436	CALL start_init_orog( iml, jml, lon_in, lat_in)
	437	ENDIF
	438	IF ( SIZE(zthe) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	439	WRITE(,)
	440	. 'STARTVAR module has been initialized to the wrong size'
	441	STOP
	442	ENDIF
	443	CALL gr_dyn_fi(1, iml, jml, nbindex,zthe, champ)
	444	CASE ('zpic')
	445	IF ( .NOT.ALLOCATED(zpic)) THEN
	446	CALL start_init_orog( iml, jml, lon_in, lat_in)
	447	ENDIF
	448	IF ( SIZE(zpic) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	449	WRITE(,)
	450	. 'STARTVAR module has been initialized to the wrong size'
	451	STOP
	452	ENDIF
	453	CALL gr_dyn_fi(1, iml, jml, nbindex,zpic, champ)
	454	CASE ('zval')
	455	IF ( .NOT.ALLOCATED(zval)) THEN
	456	CALL start_init_orog( iml, jml, lon_in, lat_in)
	457	ENDIF
	458	IF ( SIZE(zval) .NE. SIZE(lon_in)*SIZE(lat_in) ) THEN
	459	WRITE(,)
	460	. 'STARTVAR module has been initialized to the wrong size'
	461	STOP
	462	ENDIF
	463	CALL gr_dyn_fi(1, iml, jml, nbindex,zval, champ)
	464	CASE ('rads')
	465	champ(:) = 0.0
	466	CASE ('snow')
	467	champ(:) = 0.0
	468	CASE ('deltat')
	469	champ(:) = 0.0
	470	CASE ('rugmer')
	471	champ(:) = 0.001
	472	CASE ('agsno')
	473	champ(:) = 50.0
	474	CASE DEFAULT
	475	WRITE(,) 'startget_phys1d'
	476	WRITE(,) 'No rule is present to extract variable ',
	477	. varname(:LEN_TRIM(varname)),' from any data set'
	478	STOP
	479	END SELECT
	480	ELSE
	481	!
	482	! If we see tsol we catch it as we may need it for a 3D interpolation
	483	!
	484	SELECTCASE(varname)
	485	CASE ('tsol')
	486	IF ( .NOT.ALLOCATED(tsol)) THEN
	487	ALLOCATE(tsol(SIZE(lon_in),SIZE(lat_in) ))
	488	ENDIF
	489	CALL gr_fi_dyn(1, iml, jml, nbindex, champ, tsol)
	490	END SELECT
	491	ENDIF
	492	END SUBROUTINE startget_phys1d
	493	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	494	!
	495	SUBROUTINE start_init_phys( iml, jml, lon_in, lat_in)
	496	!
	497	INTEGER, INTENT(in) :: iml, jml
	498	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	499	!
	500	! LOCAL
	501	!
	502	REAL :: lev(1), date, dt
	503	INTEGER :: itau(1)
	504	INTEGER :: llm_tmp, ttm_tmp
	505	INTEGER :: i, j
	506	!
	507	CHARACTER*120 :: physfname
	508	LOGICAL :: check=.TRUE.
	509	!
	510	REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
	511	REAL, ALLOCATABLE :: var_ana(:,:), tmp_var(:,:)
	512	!
	513	physfname = 'ECPHY.nc'
	514	!
	515	IF ( check ) WRITE(,) 'Opening the surface analysis'
	516	!
	517	CALL flininfo(physfname, iml_phys, jml_phys, llm_tmp,
	518	. ttm_tmp, fid_phys)
	519	!
	520	!
	521	ALLOCATE (lat_phys(iml_phys,jml_phys))
	522	ALLOCATE (lon_phys(iml_phys,jml_phys))
	523	!
	524	CALL flinopen(physfname, .FALSE., iml_phys, jml_phys,
	525	. llm_tmp, lon_phys, lat_phys, lev, ttm_tmp,
	526	. itau, date, dt, fid_phys)
	527	!
	528	! Allocate the space we will need to get the data out of this file
	529	!
	530	ALLOCATE(var_ana(iml_phys, jml_phys))
	531	!
	532	! In case we have a file which is in degrees we do the transformation
	533	!
	534	ALLOCATE(lon_rad(iml_phys))
	535	IF ( MAXVAL(lon_phys(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	536	lon_rad(:) = lon_phys(:,1) * 2.0 * ASIN(1.0) / 180.0
	537	ELSE
	538	lon_rad(:) = lon_phys(:,1)
	539	ENDIF
	540	ALLOCATE(lat_rad(jml_phys))
	541	IF ( MAXVAL(lat_phys(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	542	lat_rad(:) = lat_phys(1,:) * 2.0 * ASIN(1.0) / 180.0
	543	ELSE
	544	lat_rad(:) = lat_phys(1,:)
	545	ENDIF
	546	!
	547	! We get the two standard varibales
	548	! Surface temperature
	549	!
	550	ALLOCATE(tsol(iml,jml))
	551	ALLOCATE(tmp_var(iml-1,jml))
	552	!
	553	!
	554	CALL flinget(fid_phys, 'ST', iml_phys, jml_phys,
	555	.llm_tmp, ttm_tmp, 1, 1, var_ana)
	556	CALL grille_m(iml_phys, jml_phys, lon_rad, lat_rad,
	557	. var_ana, iml-1, jml, lon_in, lat_in, tmp_var)
	558	CALL gr_int_dyn(tmp_var, tsol, iml-1, jml)
	559	!
	560	! Soil moisture
	561	!
	562	ALLOCATE(qsol(iml,jml))
	563	CALL flinget(fid_phys, 'CDSW', iml_phys, jml_phys,
	564	. llm_tmp, ttm_tmp, 1, 1, var_ana)
	565	CALL grille_m(iml_phys, jml_phys, lon_rad, lat_rad,
	566	. var_ana, iml-1, jml, lon_in, lat_in, tmp_var)
	567	CALL gr_int_dyn(tmp_var, qsol, iml-1, jml)
	568	!
	569	CALL flinclo(fid_phys)
	570	!
	571	END SUBROUTINE start_init_phys
	572	!
	573	!
	574	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	575	!
	576	!
	577	SUBROUTINE startget_dyn(varname, iml, jml, lon_in, lat_in,
	578	. lml, pls, workvar, champ, val_exp)
	579	!
	580	! ARGUMENTS
	581	!
	582	CHARACTER(), INTENT(in) :: varname
	583	INTEGER, INTENT(in) :: iml, jml, lml
	584	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	585	REAL, INTENT(in) :: pls(iml, jml, lml)
	586	REAL, INTENT(in) :: workvar(iml, jml, lml)
	587	REAL, INTENT(inout) :: champ(iml, jml, lml)
	588	REAL, INTENT(in) :: val_exp
	589	!
	590	! LOCAL
	591	!
	592	INTEGER :: il, ij, ii
	593	REAL :: xppn, xpps
	594	!
	595	! C'est vraiment une galere de devoir rajouter tant de commons just pour avoir les aires.
	596	! Il faudrait mettre une structure plus flexible et moins dangereuse.
	597	!
	598	#include "dimensions.h"
	599	#include "paramet.h"
	600	#include "comgeom2.h"
	601	#include "comconst.h"
	602	!
	603	! This routine only works if the variable does not exist or is constant
	604	!
	605	IF ( MINVAL(champ(:,:,:)).EQ.MAXVAL(champ(:,:,:)) .AND.
	606	. MINVAL(champ(:,:,:)).EQ.val_exp ) THEN
	607	!
	608	SELECTCASE(varname)
	609	CASE ('u')
	610	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	611	CALL start_init_dyn( iml, jml, lon_in, lat_in)
	612	ENDIF
	613	CALL start_inter_3d('U', iml, jml, lml, lon_in,
	614	. lat_in, pls, champ)
	615	DO il=1,lml
	616	DO ij=1,jml
	617	DO ii=1,iml-1
	618	champ(ii,ij,il) = champ(ii,ij,il) * cu(ii,ij)
	619	ENDDO
	620	champ(iml,ij, il) = champ(1,ij, il)
	621	ENDDO
	622	ENDDO
	623	CASE ('v')
	624	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	625	CALL start_init_dyn( iml, jml, lon_in, lat_in)
	626	ENDIF
	627	CALL start_inter_3d('V', iml, jml, lml, lon_in,
	628	. lat_in, pls, champ)
	629	DO il=1,lml
	630	DO ij=1,jml
	631	DO ii=1,iml-1
	632	champ(ii,ij,il) = champ(ii,ij,il) * cv(ii,ij)
	633	ENDDO
	634	champ(iml,ij, il) = champ(1,ij, il)
	635	ENDDO
	636	ENDDO
	637	CASE ('t')
	638	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	639	CALL start_init_dyn( iml, jml, lon_in, lat_in)
	640	ENDIF
	641	CALL start_inter_3d('TEMP', iml, jml, lml, lon_in,
	642	. lat_in, pls, champ)
	643
	644	CASE ('tpot')
	645	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	646	CALL start_init_dyn( iml, jml, lon_in, lat_in)
	647	ENDIF
	648	CALL start_inter_3d('TEMP', iml, jml, lml, lon_in,
	649	. lat_in, pls, champ)
	650	IF ( MINVAL(workvar(:,:,:)) .NE. MAXVAL(workvar(:,:,:)) )
	651	. THEN
	652	DO il=1,lml
	653	DO ij=1,jml
	654	DO ii=1,iml-1
	655	champ(ii,ij,il) = champ(ii,ij,il) * cpp
	656	. / workvar(ii,ij,il)
	657	ENDDO
	658	champ(iml,ij,il) = champ(1,ij,il)
	659	ENDDO
	660	ENDDO
	661	DO il=1,lml
	662	xppn = SUM(aire(:,1)*champ(:,1,il))/apoln
[19]	663	xpps = SUM(aire(:,jml)*champ(:,jml,il))/apols
[2]	664	champ(:,1,il) = xppn
	665	champ(:,jml,il) = xpps
	666	ENDDO
	667	ELSE
	668	WRITE(,)'Could not compute potential temperature as the'
	669	WRITE(,)'Exner function is missing or constant.'
	670	STOP
	671	ENDIF
	672	CASE ('q')
	673	IF ( .NOT.ALLOCATED(psol_dyn)) THEN
	674	CALL start_init_dyn( iml, jml, lon_in, lat_in)
	675	ENDIF
	676	CALL start_inter_3d('R', iml, jml, lml, lon_in, lat_in,
	677	. pls, champ)
	678	IF ( MINVAL(workvar(:,:,:)) .NE. MAXVAL(workvar(:,:,:)) )
	679	. THEN
	680	DO il=1,lml
	681	DO ij=1,jml
	682	DO ii=1,iml-1
	683	champ(ii,ij,il) = 0.01 * champ(ii,ij,il) *
	684	. workvar(ii,ij,il)
	685	ENDDO
	686	champ(iml,ij,il) = champ(1,ij,il)
	687	ENDDO
	688	ENDDO
	689	WHERE ( champ .LT. 0.) champ = 1.0E-10
	690	DO il=1,lml
	691	xppn = SUM(aire(:,1)*champ(:,1,il))/apoln
[19]	692	xpps = SUM(aire(:,jml)*champ(:,jml,il))/apols
[2]	693	champ(:,1,il) = xppn
	694	champ(:,jml,il) = xpps
	695	ENDDO
	696	ELSE
	697	WRITE(,)'Could not compute specific humidity as the'
	698	WRITE(,)'saturated humidity is missing or constant.'
	699	STOP
	700	ENDIF
	701	CASE DEFAULT
	702	WRITE(,) 'startget_dyn'
	703	WRITE(,) 'No rule is present to extract variable ',
	704	. varname(:LEN_TRIM(varname)),' from any data set'
	705	STOP
	706	END SELECT
	707	ENDIF
	708	END SUBROUTINE startget_dyn
	709	!
	710	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	711	!
	712	SUBROUTINE start_init_dyn( iml, jml, lon_in, lat_in)
	713	!
	714	INTEGER, INTENT(in) :: iml, jml
	715	REAL, INTENT(in) :: lon_in(iml), lat_in(jml)
	716	!
	717	! LOCAL
	718	!
	719	REAL :: lev(1), date, dt
	720	INTEGER :: itau(1)
	721	INTEGER :: i, j
	722	integer :: iret
	723	!
	724	CHARACTER*120 :: physfname
	725	LOGICAL :: check=.TRUE.
	726	!
	727	REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
	728	REAL, ALLOCATABLE :: var_ana(:,:), tmp_var(:,:), z(:,:)
	729	REAL, ALLOCATABLE :: xppn(:), xpps(:)
	730	LOGICAL :: allo
	731	!
	732	! Ce n'est pas tres pratique d'avoir a charger 3 include pour avoir la grille du modele
	733	!
	734	#include "dimensions.h"
	735	#include "paramet.h"
	736	#include "comgeom2.h"
	737	!
	738	physfname = 'ECDYN.nc'
	739	!
	740	IF ( check ) WRITE(,) 'Opening the surface analysis'
	741	!
	742	CALL flininfo(physfname, iml_dyn, jml_dyn, llm_dyn,
	743	. ttm_dyn, fid_dyn)
	744	IF ( check ) WRITE(,) 'Values read: ', iml_dyn, jml_dyn,
	745	. llm_dyn, ttm_dyn
	746	!
	747	ALLOCATE (lat_dyn(iml_dyn,jml_dyn), stat=iret)
	748	ALLOCATE (lon_dyn(iml_dyn,jml_dyn), stat=iret)
	749	ALLOCATE (lev_dyn(llm_dyn), stat=iret)
	750	!
	751	CALL flinopen(physfname, .FALSE., iml_dyn, jml_dyn, llm_dyn,
	752	. lon_dyn, lat_dyn, lev_dyn, ttm_dyn,
	753	. itau, date, dt, fid_dyn)
	754	!
	755
	756	allo = allocated (var_ana)
	757	if (allo) then
	758	DEALLOCATE(var_ana, stat=iret)
	759	endif
	760	ALLOCATE(var_ana(iml_dyn, jml_dyn), stat=iret)
	761
	762	allo = allocated (lon_rad)
	763	if (allo) then
	764	DEALLOCATE(lon_rad, stat=iret)
	765	endif
	766	ALLOCATE(lon_rad(iml_dyn), stat=iret)
	767
	768	IF ( MAXVAL(lon_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	769	lon_rad(:) = lon_dyn(:,1) * 2.0 * ASIN(1.0) / 180.0
	770	ELSE
	771	lon_rad(:) = lon_dyn(:,1)
	772	ENDIF
	773	ALLOCATE(lat_rad(jml_dyn))
	774	IF ( MAXVAL(lat_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	775	lat_rad(:) = lat_dyn(1,:) * 2.0 * ASIN(1.0) / 180.0
	776	ELSE
	777	lat_rad(:) = lat_dyn(1,:)
	778	ENDIF
	779	!
	780	ALLOCATE(z(iml, jml))
	781	ALLOCATE(tmp_var(iml-1,jml))
	782	!
	783	CALL flinget(fid_dyn, 'Z', iml_dyn, jml_dyn, 0, ttm_dyn,
	784	. 1, 1, var_ana)
	785	CALL grille_m(iml_dyn, jml_dyn , lon_rad, lat_rad, var_ana,
	786	. iml-1, jml, lon_in, lat_in, tmp_var)
	787	CALL gr_int_dyn(tmp_var, z, iml-1, jml)
	788	!
	789	ALLOCATE(psol_dyn(iml, jml))
	790	!
	791	CALL flinget(fid_dyn, 'SP', iml_dyn, jml_dyn, 0, ttm_dyn,
	792	. 1, 1, var_ana)
	793	CALL grille_m(iml_dyn, jml_dyn , lon_rad, lat_rad, var_ana,
	794	. iml-1, jml, lon_in, lat_in, tmp_var)
	795	CALL gr_int_dyn(tmp_var, psol_dyn, iml-1, jml)
	796	!
	797	IF ( .NOT.ALLOCATED(tsol)) THEN
	798	! These variables may have been allocated by the need to
	799	! create a start field for them or by the varibale
	800	! coming out of the restart file. In case we dor have it we will initialize it.
	801	!
	802	CALL start_init_phys( iml, jml, lon_in, lat_in)
	803	ELSE
	804	IF ( SIZE(tsol) .NE. SIZE(psol_dyn) ) THEN
	805	WRITE(,) 'start_init_dyn :'
	806	WRITE(,) 'The temperature field we have does not ',
	807	. 'have the right size'
	808	STOP
	809	ENDIF
	810	ENDIF
	811	IF ( .NOT.ALLOCATED(phis)) THEN
	812	!
	813	! These variables may have been allocated by the need to create a start field for them or by the varibale
	814	! coming out of the restart file. In case we dor have it we will initialize it.
	815	!
	816	CALL start_init_orog( iml, jml, lon_in, lat_in)
	817	!
	818	ELSE
	819	!
	820	IF (SIZE(phis) .NE. SIZE(psol_dyn)) THEN
	821	!
	822	WRITE(,) 'start_init_dyn :'
	823	WRITE(,) 'The orography field we have does not ',
	824	. ' have the right size'
	825	STOP
	826	ENDIF
	827	!
	828	ENDIF
	829	!
	830	! PSOL is computed in Pascals
	831	!
	832	!
	833	DO j = 1, jml
	834	DO i = 1, iml-1
	835	psol_dyn(i,j) = psol_dyn(i,j)*(1.0+(z(i,j)-phis(i,j))
	836	. /287.0/tsol(i,j))
	837	ENDDO
	838	psol_dyn(iml,j) = psol_dyn(1,j)
	839	ENDDO
	840	!
	841	!
	842	ALLOCATE(xppn(iml-1))
	843	ALLOCATE(xpps(iml-1))
	844	!
	845	DO i = 1, iml-1
	846	xppn(i) = aire( i,1) * psol_dyn( i,1)
	847	xpps(i) = aire( i,jml) * psol_dyn( i,jml)
	848	ENDDO
	849	!
	850	DO i = 1, iml
	851	psol_dyn(i,1 ) = SUM(xppn)/apoln
	852	psol_dyn(i,jml) = SUM(xpps)/apols
	853	ENDDO
	854	!
	855	RETURN
	856	!
	857	END SUBROUTINE start_init_dyn
	858	!
	859	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	860	!
	861	SUBROUTINE start_inter_3d(varname, iml, jml, lml, lon_in,
	862	. lat_in, pls_in, var3d)
	863	!
	864	! This subroutine gets a variables from a 3D file and does the interpolations needed
	865	!
	866	!
	867	! ARGUMENTS
	868	!
	869	CHARACTER() :: varname
	870	INTEGER :: iml, jml, lml
	871	REAL :: lon_in(iml), lat_in(jml), pls_in(iml, jml, lml)
	872	REAL :: var3d(iml, jml, lml)
	873	!
	874	! LOCAL
	875	!
	876	INTEGER :: ii, ij, il
	877	REAL :: bx, by
	878	REAL, ALLOCATABLE :: lon_rad(:), lat_rad(:)
	879	REAL, ALLOCATABLE :: var_tmp2d(:,:), var_tmp3d(:,:,:)
	880	REAL, ALLOCATABLE :: ax(:), ay(:), yder(:)
	881	INTEGER, ALLOCATABLE :: lind(:)
	882	!
	883	LOGICAL :: check = .TRUE.
	884	!
	885	IF ( .NOT. ALLOCATED(var_ana3d)) THEN
	886	ALLOCATE(var_ana3d(iml_dyn, jml_dyn, llm_dyn))
	887	ENDIF
	888	!
	889	!
	890	IF ( check) WRITE(,) 'Going into flinget to extract the 3D ',
	891	. ' field.', fid_dyn
	892	IF ( check) WRITE(,) fid_dyn, varname, iml_dyn, jml_dyn,
	893	. llm_dyn,ttm_dyn
	894	!
	895	CALL flinget(fid_dyn, varname, iml_dyn, jml_dyn, llm_dyn,
	896	. ttm_dyn, 1, 1, var_ana3d)
	897	!
	898	IF ( check) WRITE(,) 'Allocating space for the interpolation',
	899	. iml, jml, llm_dyn
	900	!
	901	ALLOCATE(lon_rad(iml_dyn))
	902	IF ( MAXVAL(lon_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	903	lon_rad(:) = lon_dyn(:,1) * 2.0 * ASIN(1.0) / 180.0
	904	ELSE
	905	lon_rad(:) = lon_dyn(:,1)
	906	ENDIF
	907	ALLOCATE(lat_rad(jml_dyn))
	908	IF ( MAXVAL(lat_dyn(:,:)) .GT. 2.0 * ASIN(1.0) ) THEN
	909	lat_rad(:) = lat_dyn(1,:) * 2.0 * ASIN(1.0) / 180.0
	910	ELSE
	911	lat_rad(:) = lat_dyn(1,:)
	912	ENDIF
	913	!
	914	ALLOCATE(var_tmp2d(iml-1, jml))
	915	ALLOCATE(var_tmp3d(iml, jml, llm_dyn))
	916	ALLOCATE(ax(llm_dyn))
	917	ALLOCATE(ay(llm_dyn))
	918	ALLOCATE(yder(llm_dyn))
	919	ALLOCATE(lind(llm_dyn))
	920	!
	921	DO il=1,llm_dyn
	922	!
	923	CALL grille_m(iml_dyn, jml_dyn, lon_rad, lat_rad,
	924	.var_ana3d(:,:,il), iml-1, jml, lon_in, lat_in, var_tmp2d)
	925	!
	926	CALL gr_int_dyn(var_tmp2d, var_tmp3d(:,:,il), iml-1, jml)
	927	!
	928	ENDDO
	929	!
	930	! IF needed we return the vertical axis. The spline interpolation
	931	! Requires the coordinate to be in increasing order.
	932	!
	933	IF ( lev_dyn(1) .LT. lev_dyn(llm_dyn)) THEN
	934	DO il=1,llm_dyn
	935	lind(il) = il
	936	ENDDO
	937	ELSE
	938	DO il=1,llm_dyn
	939	lind(il) = llm_dyn-il+1
	940	ENDDO
	941	ENDIF
	942	!
	943	DO ij=1,jml
	944	DO ii=1,iml-1
	945	!
	946	ax(:) = lev_dyn(lind(:)) * 100
	947	ay(:) = var_tmp3d(ii, ij, lind(:))
	948	!
	949	CALL SPLINE(ax, ay, llm_dyn, 1.e30, 1.e30, yder)
	950	!
	951	DO il=1,lml
	952	bx = pls_in(ii, ij, il)
	953	CALL SPLINT(ax, ay, yder, llm_dyn, bx, by)
	954	var3d(ii, ij, il) = by
	955	ENDDO
	956	!
	957	ENDDO
	958	var3d(iml, ij, :) = var3d(1, ij, :)
	959	ENDDO
	960
	961	DEALLOCATE(lon_rad)
	962	DEALLOCATE(lat_rad)
	963	DEALLOCATE(var_tmp2d)
	964	DEALLOCATE(var_tmp3d)
	965	DEALLOCATE(ax)
	966	DEALLOCATE(ay)
	967	DEALLOCATE(yder)
	968	DEALLOCATE(lind)
	969
	970	!
	971	RETURN
	972	!
	973	END SUBROUTINE start_inter_3d
	974	!
	975	END MODULE startvar

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

Download in other formats:

Original Format