Context Navigation

datareadnc.F @ 161

Last change on this file since 161 was 146, checked in by emillour, 14 years ago

Mars GCM:

minor bug fix in lect_start_archive.F (using wrong surface temp. array).

+ swiched output messages to english and added that tracers not found
in file must be initialized by user.

minor bug fix in datareadnc.F : 'datafile' path must be initialized.

File size: 9.4 KB

Line
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	! to use 'getin'
45	use ioipsl_getincom
46
47	implicit none
48
49	#include "dimensions.h"
50	#include "paramet.h"
51	#include "comgeom.h"
52	#include "comconst.h"
53	#include "netcdf.inc"
54	#include "datafile.h"
55
56	c=======================================================================
57	c Declarations:
58	C=======================================================================
59
60	INTEGER imd,jmd,imdp1,jmdp1
61	parameter (imd=360,jmd=179,imdp1=361,jmdp1=180)
62
63	INTEGER iimp1
64	parameter (iimp1=iim+1-1/iim)
65
66	CHARACTER relief*3
67
68	REAL zdata(imd*jmdp1)
69	REAL zdataS(imdp1*jmdp1)
70	REAL pfield(iimp1*jjp1)
71
72	REAL alb(iimp1*jjp1)
73	REAL ith(iimp1*jjp1)
74	REAL phisinit(iimp1*jjp1)
75
76	REAL zmea(imdp1*jmdp1)
77	REAL zstd(imdp1*jmdp1)
78	REAL zsig(imdp1*jmdp1)
79	REAL zgam(imdp1*jmdp1)
80	REAL zthe(imdp1*jmdp1)
81
82	INTEGER ierr
83
84	INTEGER unit,nvarid
85
86	INTEGER i,j,k
87
88	INTEGER klatdat,ngridmxgdat
89	PARAMETER (klatdat=180,ngridmxgdat=360)
90
91	c on passe une grille en rlonu rlatv et im+1 jm a interp_horiz)
92
93	REAL longitude(imd),latitude(jmdp1) ! Pour lecture des donnees
94	REAL rlonud(imdp1),rlatvd(jmd)
95
96	CHARACTER*20 string
97	DIMENSION string(4)
98
99
100	#include "lmdstd.h"
101	#include "fxyprim.h"
102
103	pi=2.*ASIN(1.)
104
105	c=======================================================================
106	c rlonud, rlatvd
107	c=======================================================================
108
109	c-----------------------------------------------------------------------
110	c Lecture NetCDF des donnees latitude et longitude
111	c-----------------------------------------------------------------------
112	write(,) 'datareadnc: opening file surface.nc'
113
114	datafile="/u/forget/WWW/datagcm/datafile" ! default path to surface.nc
115	call getin("datadir",datafile) ! but users may specify another path
116
117	ierr = NF_OPEN (trim(datafile)//'/surface.nc',
118	& NF_NOWRITE,unit)
119	IF (ierr.NE.NF_NOERR) THEN
120	write(,)'Error : cannot open file surface.nc '
121	write(,)'(in phymars/datareadnc.F)'
122	write(,)'It should be in :',trim(datafile),'/'
123	write(,)'1) You can set this path in the callphys.def file:'
124	write(,)' datadir=/path/to/the/datafiles'
125	write(,)'2) If necessary surface.nc (and other datafiles)'
126	write(,)' can be obtained online on:'
127	write(,)' http://www.lmd.jussieu.fr/~forget/datagcm/datafile'
128	CALL ABORT
129	ENDIF
130
131	c
132	c Lecture des latitudes (coordonnees):
133	c
134	ierr = NF_INQ_VARID (unit, "latitude", nvarid)
135	#ifdef NC_DOUBLE
136	ierr = NF_GET_VAR_DOUBLE(unit, nvarid, latitude)
137	#else
138	ierr = NF_GET_VAR_REAL(unit, nvarid, latitude)
139	#endif
140	c
141	c Lecture des longitudes (coordonnees):
142	c
143	ierr = NF_INQ_VARID (unit, "longitude", nvarid)
144	#ifdef NC_DOUBLE
145	ierr = NF_GET_VAR_DOUBLE(unit, nvarid, longitude)
146	#else
147	ierr = NF_GET_VAR_REAL(unit, nvarid, longitude)
148	#endif
149
150	c-----------------------------------------------------------------------
151	c Passage au format boites scalaires
152	c-----------------------------------------------------------------------
153
154	c-----------------------------------------------------------------------
155	c longitude(imd) --> rlonud(imdp1)
156	c-----------------------------------------------------------------------
157
158	c Passage en coordonnees boites scalaires et en radian
159	do i=1,imd
160	rlonud(i)=(longitude(i)+.5)*pi/180.
161	enddo
162
163	c Repetition de la valeur im+1
164	rlonud(imdp1)=rlonud(1) + 2*pi
165
166	c-----------------------------------------------------------------------
167	c latitude(jmdp1) --> rlonvd(jmd)
168	c-----------------------------------------------------------------------
169
170	c Passage en coordonnees boites scalaires et en radian
171	do j=1,jmd
172	rlatvd(j)=(latitude(j)-.5)*pi/180.
173	enddo
174
175	c=======================================================================
176	c lecture NetCDF de albedo, thermal, relief, zdtm (pour francois Lott)
177	c=======================================================================
178
179	string(1) = 'albedo'
180	string(2) = 'thermal'
181	if (relief.ne.'pla') then
182	write(,) ' MOLA topography'
183	relief = 'MOL'
184	string(3) = 'z'//relief
185	else
186	string(3) = 'zMOL' ! pour qu''il lise qqchose sur le fichier
187	! remise a 0 derriere
188	endif
189	string(4) = 'zMOL' ! lecture pour calcul topog. sous-maille
190
191
192	DO k=1,4
193	write(,) 'string',k,string(k)
194
195	c-----------------------------------------------------------------------
196	c initialisation
197	c-----------------------------------------------------------------------
198	call initial0(iimp1*jjp1,pfield)
199	call initial0(imd*jmdp1,zdata)
200	call initial0(imdp1*jmdp1,zdataS)
201
202	c-----------------------------------------------------------------------
203	c Lecture NetCDF
204	c-----------------------------------------------------------------------
205
206	ierr = NF_INQ_VARID (unit, string(k), nvarid)
207	#ifdef NC_DOUBLE
208	ierr = NF_GET_VAR_DOUBLE(unit, nvarid, zdata)
209	#else
210	ierr = NF_GET_VAR_REAL(unit, nvarid, zdata)
211	#endif
212
213	c-----------------------------------------------------------------------
214	c Cas particulier "Francois Lott" ( k=4 ) (relief sous-maille)
215	c-----------------------------------------------------------------------
216	if (k.eq.4) then
217
218	call multscal(imd*jmdp1,zdata,1000.,zdata)
219	call multscal(imd,longitude,pi/180.,longitude)
220	call multscal(jmdp1,latitude,pi/180.,latitude)
221
222	call grid_noro1(360, 180, longitude, latitude, zdata,
223	. iim, jjp1, rlonv, rlatu, zmea,zstd,zsig,zgam,zthe)
224
225	CALL dump2d(iip1,jjp1,zmea,'zmea')
226	CALL dump2d(iip1,jjp1,zstd,'zstd')
227	CALL dump2d(iip1,jjp1,zsig,'zsig')
228	CALL dump2d(iip1,jjp1,zgam,'zgam')
229	CALL dump2d(iip1,jjp1,zthe,'zthe')
230
231	endif
232
233	c-----------------------------------------------------------------------
234	c Passage de zdata en grille (imdp1 jmdp1)
235	c-----------------------------------------------------------------------
236	do j=1,jmdp1
237	do i=1,imd
238	zdataS(i+imdp1(j-1)) = zdata(i+ngridmxgdat(j-1))
239	enddo
240	zdataS(imdp1+imdp1(j-1)) = zdata(1+ngridmxgdat(j-1))
241	enddo
242
243	c-----------------------------------------------------------------------
244	c Interpolation
245	c-----------------------------------------------------------------------
246	call interp_horiz(zdataS,pfield,imd,jmd,
247	. iim, jjm,1,rlonud,rlatvd,rlonu,rlatv)
248
249	c-----------------------------------------------------------------------
250	c Periodicite
251	c-----------------------------------------------------------------------
252
253	do j=1,jjp1
254	pfield(iimp1j) = pfield(1+iimp1(j-1))
255	enddo
256
257	c-----------------------------------------------------------------------
258	c Sauvegarde des champs
259	c-----------------------------------------------------------------------
260
261	if (k.eq.1) then ! albedo
262	do i=1,iimp1*jjp1
263	alb(i) = pfield(i)
264	enddo
265	elseif (k.eq.2) then ! thermal
266	do i=1,iimp1*jjp1
267	ith(i) = pfield(i)
268	enddo
269	elseif (k.eq.3) then ! relief
270	if (relief.eq.'pla') then
271	call initial0(iimp1*jjp1,phisinit)
272	else
273	do i=1,iimp1*jjp1
274	phisinit(i) = pfield(i)
275	enddo
276	endif
277	endif
278
279	ENDDO
280
281	c-----------------------------------------------------------------------
282	c Traitement Phisinit
283	c-----------------------------------------------------------------------
284
285	DO i=1,iimp1*jjp1
286	phisinit(i)=1000.*phisinit(i)
287	ENDDO
288	CALL dump2d(iimp1,jjp1,phisinit,'Altitude in m')
289	CALL multscal(iimp1*jjp1,phisinit,g,phisinit)
290
291	c-----------------------------------------------------------------------
292	c FIN
293	c-----------------------------------------------------------------------
294
295	END

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

Download in other formats:

Original Format