Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

datareadnc.F @ 57

Last change on this file since 57 was 57, checked in by aslmd, 14 years ago
mineur LMD_MM_MARS: ajout du GCM ancienne physique, systeme maintenant complet sur SVN (ne manque que la base de donnees d'etats initiaux)
File size: 9.2 KB

Rev	Line
[57]	1	c=======================================================================
	2	SUBROUTINE datareadnc(relief,phisinit,alb,ith,
	3	. zmea,zstd,zsig,zgam,zthe)
	4	c=======================================================================
	5	c
	6	c
	7	c Author: F. Hourdin 01/1997
	8	c -------
	9	c
	10	c Object: To read data from Martian surface to use in a GCM
	11	c ------ from NetCDF file "surface.nc"
	12	c
	13	c
	14	c Arguments:
	15	c ----------
	16	c
	17	c Inputs:
	18	c ------
	19	c
	20	c Outputs:
	21	c --------
	22	c
	23	c=======================================================================
	24	c donnees ALBEDO, INERTIE THERMIQUE, RELIEF:
	25	c
	26	c Ces donnees sont au format NetCDF dans le fichier "surface.nc"
	27	c
	28	c 360 valeurs en longitude (de -179.5 a 179.5)
	29	c 180 valeurs en latitudes (de 89.5 a -89.5)
	30	c
	31	c Pour les passer au format de la grille, on utilise "interp_horiz.F"
	32	c
	33	c Il faut donc que ces donnees soient au format grille scalaire
	34	c (imold+1 jmold+1)
	35	c avec passage des coordonnees de la "boite" (rlonu, rlatv)
	36	c
	37	c On prend imd (d pour donnees!)
	38	c imd = 360 avec copie de la 1ere valeur sur la imd+1
	39	c (rlonud de -179 a -181)
	40	c jmd = 179
	41	c (rlatvd de 89 a -89)
	42	c=======================================================================
	43
	44	implicit none
	45
	46	#include "dimensions.h"
	47	#include "paramet.h"
	48	#include "comgeom.h"
	49	#include "comconst.h"
	50	#include "netcdf.inc"
	51	#include "datafile.h"
	52
	53	c=======================================================================
	54	c Declarations:
	55	C=======================================================================
	56
	57	INTEGER imd,jmd,imdp1,jmdp1
	58	parameter (imd=360,jmd=179,imdp1=361,jmdp1=180)
	59
	60	INTEGER iimp1
	61	parameter (iimp1=iim+1-1/iim)
	62
	63	CHARACTER relief*3
	64
	65	REAL zdata(imd*jmdp1)
	66	REAL zdataS(imdp1*jmdp1)
	67	REAL pfield(iimp1*jjp1)
	68
	69	REAL alb(iimp1*jjp1)
	70	REAL ith(iimp1*jjp1)
	71	REAL phisinit(iimp1*jjp1)
	72
	73	REAL zmea(imdp1*jmdp1)
	74	REAL zstd(imdp1*jmdp1)
	75	REAL zsig(imdp1*jmdp1)
	76	REAL zgam(imdp1*jmdp1)
	77	REAL zthe(imdp1*jmdp1)
	78
	79	INTEGER lnblnk, ierr
	80	EXTERNAL lnblnk
	81
	82	INTEGER unit,nvarid
	83
	84	INTEGER i,j,k
	85
	86	INTEGER klatdat,ngridmxgdat
	87	PARAMETER (klatdat=180,ngridmxgdat=360)
	88
	89	c on passe une grille en rlonu rlatv et im+1 jm a interp_horiz)
	90
	91	REAL longitude(imd),latitude(jmdp1) ! Pour lecture des donnees
	92	REAL rlonud(imdp1),rlatvd(jmd)
	93
	94	CHARACTER*20 string
	95	DIMENSION string(4)
	96
	97	CHARACTER*80 file
	98
	99	#include "lmdstd.h"
	100	#include "fxyprim.h"
	101
	102	pi=2.*ASIN(1.)
	103
	104	c=======================================================================
	105	c rlonud, rlatvd
	106	c=======================================================================
	107
	108	c-----------------------------------------------------------------------
	109	c Lecture NetCDF des donnees latitude et longitude
	110	c-----------------------------------------------------------------------
	111	write(,) 'ouverture du fichier surface.nc'
	112
	113	ierr = NF_OPEN (datafile(1:lnblnk(datafile))//'/surface.nc',
	114	& NF_NOWRITE,unit)
	115	IF (ierr.NE.NF_NOERR) THEN
	116	write(,)'Error : cannot open file surface.nc '
	117	write(,)'(in phymars/datareadnc.F)'
	118	write(,)'It should be in :',datafile(1:lnblnk(datafile)),'/'
	119	write(,)'1) You can change this directory address in '
	120	write(,)' file phymars/datafile.h'
	121	write(,)'2) If necessary surface.nc (and other datafiles)'
	122	write(,)' can be obtained online on:'
	123	write(,)' http://www.lmd.jussieu.fr/~forget/datagcm/datafile'
	124	CALL ABORT
	125	ENDIF
	126
	127	c
	128	c Lecture des latitudes (coordonnees):
	129	c
	130	ierr = NF_INQ_VARID (unit, "latitude", nvarid)
	131	#ifdef NC_DOUBLE
	132	ierr = NF_GET_VAR_DOUBLE(unit, nvarid, latitude)
	133	#else
	134	ierr = NF_GET_VAR_REAL(unit, nvarid, latitude)
	135	#endif
	136	c
	137	c Lecture des longitudes (coordonnees):
	138	c
	139	ierr = NF_INQ_VARID (unit, "longitude", nvarid)
	140	#ifdef NC_DOUBLE
	141	ierr = NF_GET_VAR_DOUBLE(unit, nvarid, longitude)
	142	#else
	143	ierr = NF_GET_VAR_REAL(unit, nvarid, longitude)
	144	#endif
	145
	146	c-----------------------------------------------------------------------
	147	c Passage au format boites scalaires
	148	c-----------------------------------------------------------------------
	149
	150	c-----------------------------------------------------------------------
	151	c longitude(imd) --> rlonud(imdp1)
	152	c-----------------------------------------------------------------------
	153
	154	c Passage en coordonnees boites scalaires et en radian
	155	do i=1,imd
	156	rlonud(i)=(longitude(i)+.5)*pi/180.
	157	enddo
	158
	159	c Repetition de la valeur im+1
	160	rlonud(imdp1)=rlonud(1) + 2*pi
	161
	162	c-----------------------------------------------------------------------
	163	c latitude(jmdp1) --> rlonvd(jmd)
	164	c-----------------------------------------------------------------------
	165
	166	c Passage en coordonnees boites scalaires et en radian
	167	do j=1,jmd
	168	rlatvd(j)=(latitude(j)-.5)*pi/180.
	169	enddo
	170
	171	c=======================================================================
	172	c lecture NetCDF de albedo, thermal, relief, zdtm (pour francois Lott)
	173	c=======================================================================
	174
	175	string(1) = 'albedo'
	176	string(2) = 'thermal'
	177	if (relief.ne.'pla') then
	178	write(,) ' La topographie est celle de MOLA'
	179	relief = 'MOL'
	180	string(3) = 'z'//relief
	181	else
	182	string(3) = 'zMOL' ! pour qu''il lise qqchose sur le fichier
	183	! remise a 0 derriere
	184	endif
	185	string(4) = 'zMOL' ! lecture pour calcul topog. sous-maille
	186
	187
	188	DO k=1,4
	189	write(,) 'string',k,string(k)
	190
	191	c-----------------------------------------------------------------------
	192	c initialisation
	193	c-----------------------------------------------------------------------
	194	call initial0(iimp1*jjp1,pfield)
	195	call initial0(imd*jmdp1,zdata)
	196	call initial0(imdp1*jmdp1,zdataS)
	197
	198	c-----------------------------------------------------------------------
	199	c Lecture NetCDF
	200	c-----------------------------------------------------------------------
	201
	202	ierr = NF_INQ_VARID (unit, string(k), nvarid)
	203	#ifdef NC_DOUBLE
	204	ierr = NF_GET_VAR_DOUBLE(unit, nvarid, zdata)
	205	#else
	206	ierr = NF_GET_VAR_REAL(unit, nvarid, zdata)
	207	#endif
	208
	209	c-----------------------------------------------------------------------
	210	c Cas particulier "Francois Lott" ( k=4 ) (relief sous-maille)
	211	c-----------------------------------------------------------------------
	212	if (k.eq.4) then
	213
	214	call multscal(imd*jmdp1,zdata,1000.,zdata)
	215	call multscal(imd,longitude,pi/180.,longitude)
	216	call multscal(jmdp1,latitude,pi/180.,latitude)
	217
	218	call grid_noro1(360, 180, longitude, latitude, zdata,
	219	. iim, jjp1, rlonv, rlatu, zmea,zstd,zsig,zgam,zthe)
	220
	221	CALL dump2d(iip1,jjp1,zmea,'zmea')
	222	CALL dump2d(iip1,jjp1,zstd,'zstd')
	223	CALL dump2d(iip1,jjp1,zsig,'zsig')
	224	CALL dump2d(iip1,jjp1,zgam,'zgam')
	225	CALL dump2d(iip1,jjp1,zthe,'zthe')
	226
	227	endif
	228
	229	c-----------------------------------------------------------------------
	230	c Passage de zdata en grille (imdp1 jmdp1)
	231	c-----------------------------------------------------------------------
	232	do j=1,jmdp1
	233	do i=1,imd
	234	zdataS(i+imdp1(j-1)) = zdata(i+ngridmxgdat(j-1))
	235	enddo
	236	zdataS(imdp1+imdp1(j-1)) = zdata(1+ngridmxgdat(j-1))
	237	enddo
	238
	239	c-----------------------------------------------------------------------
	240	c Interpolation
	241	c-----------------------------------------------------------------------
	242	call interp_horiz(zdataS,pfield,imd,jmd,
	243	. iim, jjm,1,rlonud,rlatvd,rlonu,rlatv)
	244
	245	c-----------------------------------------------------------------------
	246	c Periodicite
	247	c-----------------------------------------------------------------------
	248
	249	do j=1,jjp1
	250	pfield(iimp1j) = pfield(1+iimp1(j-1))
	251	enddo
	252
	253	c-----------------------------------------------------------------------
	254	c Sauvegarde des champs
	255	c-----------------------------------------------------------------------
	256
	257	if (k.eq.1) then ! albedo
	258	do i=1,iimp1*jjp1
	259	alb(i) = pfield(i)
	260	enddo
	261	elseif (k.eq.2) then ! thermal
	262	do i=1,iimp1*jjp1
	263	ith(i) = pfield(i)
	264	enddo
	265	elseif (k.eq.3) then ! relief
	266	if (relief.eq.'pla') then
	267	call initial0(iimp1*jjp1,phisinit)
	268	else
	269	do i=1,iimp1*jjp1
	270	phisinit(i) = pfield(i)
	271	enddo
	272	endif
	273	endif
	274
	275	ENDDO
	276
	277	c-----------------------------------------------------------------------
	278	c Traitement Phisinit
	279	c-----------------------------------------------------------------------
	280
	281	DO i=1,iimp1*jjp1
	282	phisinit(i)=1000.*phisinit(i)
	283	ENDDO
	284	CALL dump2d(iimp1,jjp1,phisinit,'Altitude en m')
	285	CALL multscal(iimp1*jjp1,phisinit,g,phisinit)
	286
	287	c-----------------------------------------------------------------------
	288	c FIN
	289	c-----------------------------------------------------------------------
	290
	291	END

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

Download in other formats:

Original Format