Context Navigation

iniwrite.F @ 1532

Last change on this file since 1532 was 1532, checked in by emillour, 9 years ago

Mars GCM:

Some fixes for buggy outputs in 1D introduced by previous code modifications.
made statto.h a module.
ecritphy in dyn3d/control_mod.F90 should be an integer, not a real.

File size: 11.7 KB

Line
1	SUBROUTINE iniwrite(nid,idayref,phis,area,nbplon,nbplat)
2
3	use comsoil_h, only: mlayer, nsoilmx
4	USE comcstfi_h, only: g, mugaz, omeg, rad, rcp, pi
5	USE comvert_mod, ONLY: ap,bp,aps,bps,pseudoalt
6	USE logic_mod, ONLY: fxyhypb,ysinus
7	USE serre_mod, ONLY: clon,clat,grossismx,grossismy,dzoomx,dzoomy
8	USE time_phylmdz_mod, ONLY: hour_ini, daysec, dtphys
9	USE ener_mod, ONLY: etot0,ptot0,ztot0,stot0,ang0
10	USE regular_lonlat_mod, ONLY: lon_reg, lat_reg
11	USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat, nbp_lev
12	IMPLICIT NONE
13
14	c=======================================================================
15	c
16	c Auteur: L. Fairhead , P. Le Van, Y. Wanherdrick, F. Forget
17	c -------
18	c
19	c Objet:
20	c ------
21	c
22	c 'Initialize' the diagfi.nc file: write down dimensions as well
23	c as time-independent fields (e.g: geopotential, mesh area, ...)
24	c
25	c=======================================================================
26	c-----------------------------------------------------------------------
27	c Declarations:
28	c -------------
29
30	include "netcdf.inc"
31
32	c Arguments:
33	c ----------
34
35	integer,intent(in) :: nid ! NetCDF file ID
36	INTEGER*4,intent(in) :: idayref ! date (initial date for this run)
37	real,intent(in) :: phis(nbplon,nbp_lat) ! surface geopotential
38	real,intent(in) :: area(nbplon,nbp_lat) ! mesh area (m2)
39	integer,intent(in) :: nbplon,nbplat ! sizes of area and phis arrays
40
41	c Local:
42	c ------
43	INTEGER length,l
44	parameter (length = 100)
45	REAL tab_cntrl(length) ! run parameters are stored in this array
46	INTEGER ierr
47	REAl,ALLOCATABLE :: lon_reg_ext(:) ! extended longitudes
48
49	integer :: nvarid,idim_index,idim_rlonv
50	integer :: idim_rlatu,idim_llmp1,idim_llm
51	integer :: idim_nsoilmx ! "subsurface_layers" dimension ID #
52	integer, dimension(2) :: id
53	c-----------------------------------------------------------------------
54
55	IF (nbp_lon*nbp_lat==1) THEN
56	! 1D model
57	ALLOCATE(lon_reg_ext(1))
58	ELSE
59	! 3D model
60	ALLOCATE(lon_reg_ext(nbp_lon+1))
61	ENDIF
62
63	DO l=1,length
64	tab_cntrl(l)=0.
65	ENDDO
66	tab_cntrl(1) = real(nbp_lon)
67	tab_cntrl(2) = real(nbp_lat-1)
68	tab_cntrl(3) = real(nbp_lev)
69	tab_cntrl(4) = real(idayref)
70	tab_cntrl(5) = rad
71	tab_cntrl(6) = omeg
72	tab_cntrl(7) = g
73	tab_cntrl(8) = mugaz
74	tab_cntrl(9) = rcp
75	tab_cntrl(10) = daysec
76	tab_cntrl(11) = dtphys
77	tab_cntrl(12) = etot0
78	tab_cntrl(13) = ptot0
79	tab_cntrl(14) = ztot0
80	tab_cntrl(15) = stot0
81	tab_cntrl(16) = ang0
82
83	tab_cntrl(27) = hour_ini
84	c
85	c .......... P.Le Van ( ajout le 8/04/96 ) .........
86	c ..... parametres pour le zoom ......
87	tab_cntrl(17) = clon
88	tab_cntrl(18) = clat
89	tab_cntrl(19) = grossismx
90	tab_cntrl(20) = grossismy
91	c
92	c ..... ajout le 6/05/97 et le 15/10/97 .......
93	c
94	IF ( fxyhypb ) THEN
95	tab_cntrl(21) = 1.
96	tab_cntrl(22) = dzoomx
97	tab_cntrl(23) = dzoomy
98	ELSE
99	tab_cntrl(21) = 0.
100	tab_cntrl(22) = dzoomx
101	tab_cntrl(23) = dzoomy
102	tab_cntrl(24) = 0.
103	IF( ysinus ) tab_cntrl(24) = 1.
104	ENDIF
105
106	c .........................................................
107
108	! Define dimensions
109
110	ierr = NF_REDEF (nid)
111
112	ierr = NF_DEF_DIM (nid, "index", length, idim_index)
113	! ierr = NF_DEF_DIM (nid, "rlonu", iip1, idim_rlonu)
114	ierr = NF_DEF_DIM (nid, "latitude", nbp_lat, idim_rlatu)
115	IF (nbp_lon*nbp_lat==1) THEN
116	ierr = NF_DEF_DIM (nid, "longitude", 1, idim_rlonv)
117	ELSE
118	ierr = NF_DEF_DIM (nid, "longitude", nbp_lon+1, idim_rlonv)
119	ENDIF
120	! ierr = NF_DEF_DIM (nid, "rlatv", jjm, idim_rlatv)
121	ierr = NF_DEF_DIM (nid, "interlayer", (nbp_lev+1), idim_llmp1)
122	ierr = NF_DEF_DIM (nid, "altitude", nbp_lev, idim_llm)
123	ierr = NF_DEF_DIM (nid,"subsurface_layers",nsoilmx,idim_nsoilmx)
124	c
125	ierr = NF_ENDDEF(nid)
126
127	c Contol parameters for this run
128	ierr = NF_REDEF (nid)
129	#ifdef NC_DOUBLE
130	ierr = NF_DEF_VAR (nid, "controle", NF_DOUBLE, 1,
131	. idim_index,nvarid)
132	#else
133	ierr = NF_DEF_VAR (nid, "controle", NF_FLOAT, 1,
134	. idim_index,nvarid)
135	#endif
136	ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 18,
137	. "Control parameters")
138	ierr = NF_ENDDEF(nid)
139	#ifdef NC_DOUBLE
140	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tab_cntrl)
141	#else
142	ierr = NF_PUT_VAR_REAL (nid,nvarid,tab_cntrl)
143	#endif
144
145	c --------------------------
146	c longitudes and latitudes
147	!
148	! ierr = NF_REDEF (nid)
149	!#ifdef NC_DOUBLE
150	! ierr = NF_DEF_VAR (nid, "rlonu", NF_DOUBLE, 1, idim_rlonu,nvarid)
151	!#else
152	! ierr = NF_DEF_VAR (nid, "rlonu", NF_FLOAT, 1, idim_rlonu,nvarid)
153	!#endif
154	! ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 21,
155	! . "Longitudes at u nodes")
156	! ierr = NF_ENDDEF(nid)
157	!#ifdef NC_DOUBLE
158	! ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlonu/pi*180)
159	!#else
160	! ierr = NF_PUT_VAR_REAL (nid,nvarid,rlonu/pi*180)
161	!#endif
162	c
163	c --------------------------
164	ierr = NF_REDEF (nid)
165	#ifdef NC_DOUBLE
166	ierr =NF_DEF_VAR(nid, "latitude", NF_DOUBLE, 1, idim_rlatu,nvarid)
167	#else
168	ierr =NF_DEF_VAR (nid, "latitude", NF_FLOAT, 1, idim_rlatu,nvarid)
169	#endif
170	ierr =NF_PUT_ATT_TEXT(nid,nvarid,'units',13,"degrees_north")
171	ierr = NF_PUT_ATT_TEXT (nid,nvarid,"long_name", 14,
172	. "North latitude")
173	ierr = NF_ENDDEF(nid)
174	#ifdef NC_DOUBLE
175	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,lat_reg/pi*180)
176	#else
177	ierr = NF_PUT_VAR_REAL (nid,nvarid,lat_reg/pi*180)
178	#endif
179	c
180	c --------------------------
181
182	lon_reg_ext(1:nbp_lon)=lon_reg(1:nbp_lon)
183	IF (nbp_lon*nbp_lat/=1) THEN
184	! In 3D, add extra redundant point (180 degrees,
185	! since lon_reg starts at -180)
186	lon_reg_ext(nbp_lon+1)=-lon_reg_ext(1)
187	ENDIF
188
189	ierr = NF_REDEF (nid)
190	#ifdef NC_DOUBLE
191	ierr =NF_DEF_VAR(nid,"longitude", NF_DOUBLE, 1, idim_rlonv,nvarid)
192	#else
193	ierr = NF_DEF_VAR(nid,"longitude", NF_FLOAT, 1, idim_rlonv,nvarid)
194	#endif
195	ierr = NF_PUT_ATT_TEXT (nid,nvarid,"long_name", 14,
196	. "East longitude")
197	ierr = NF_PUT_ATT_TEXT(nid,nvarid,'units',12,"degrees_east")
198	ierr = NF_ENDDEF(nid)
199	#ifdef NC_DOUBLE
200	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,lon_reg_ext/pi*180)
201	#else
202	ierr = NF_PUT_VAR_REAL (nid,nvarid,lon_reg_ext/pi*180)
203	#endif
204	c
205	c --------------------------
206	ierr = NF_REDEF (nid)
207	#ifdef NC_DOUBLE
208	ierr = NF_DEF_VAR (nid, "altitude", NF_DOUBLE, 1,
209	. idim_llm,nvarid)
210	#else
211	ierr = NF_DEF_VAR (nid, "altitude", NF_FLOAT, 1,
212	. idim_llm,nvarid)
213	#endif
214	ierr = NF_PUT_ATT_TEXT (nid,nvarid,"long_name",10,"pseudo-alt")
215	ierr = NF_PUT_ATT_TEXT (nid,nvarid,'units',2,"km")
216	ierr = NF_PUT_ATT_TEXT (nid,nvarid,'positive',2,"up")
217
218	ierr = NF_ENDDEF(nid)
219	#ifdef NC_DOUBLE
220	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,pseudoalt)
221	#else
222	ierr = NF_PUT_VAR_REAL (nid,nvarid,pseudoalt)
223	#endif
224	c
225	c --------------------------
226	! ierr = NF_REDEF (nid)
227	!#ifdef NC_DOUBLE
228	! ierr = NF_DEF_VAR (nid, "rlatv", NF_DOUBLE, 1, idim_rlatv,nvarid)
229	!#else
230	! ierr = NF_DEF_VAR (nid, "rlatv", NF_FLOAT, 1, idim_rlatv,nvarid)
231	!#endif
232	! ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 20,
233	! . "Latitudes at v nodes")
234	! ierr = NF_ENDDEF(nid)
235	!#ifdef NC_DOUBLE
236	! ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlatv/pi*180)
237	!#else
238	! ierr = NF_PUT_VAR_REAL (nid,nvarid,rlatv/pi*180)
239	!#endif
240	c
241	c --------------------------
242	c Vertical levels
243	call def_var(nid,"aps","hybrid pressure at midlayers ","Pa",
244	. 1,idim_llm,nvarid,ierr)
245	#ifdef NC_DOUBLE
246	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,aps)
247	#else
248	ierr = NF_PUT_VAR_REAL (nid,nvarid,aps)
249	#endif
250
251	call def_var(nid,"bps","hybrid sigma at midlayers"," ",
252	. 1,idim_llm,nvarid,ierr)
253	#ifdef NC_DOUBLE
254	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,bps)
255	#else
256	ierr = NF_PUT_VAR_REAL (nid,nvarid,bps)
257	#endif
258
259	call def_var(nid,"ap","hybrid pressure at interlayers","Pa",
260	. 1,idim_llmp1,nvarid,ierr)
261	#ifdef NC_DOUBLE
262	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,ap)
263	#else
264	ierr = NF_PUT_VAR_REAL (nid,nvarid,ap)
265	#endif
266
267	call def_var(nid,"bp","hybrid sigma at interlayers"," ",
268	. 1,idim_llmp1,nvarid,ierr)
269	#ifdef NC_DOUBLE
270	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,bp)
271	#else
272	ierr = NF_PUT_VAR_REAL (nid,nvarid,bp)
273	#endif
274
275	!-------------------------------
276	! (soil) depth variable mlayer() (known from comsoil_h)
277	!-------------------------------
278	ierr=NF_REDEF (nid) ! Enter NetCDF (re-)define mode
279	! define variable
280	#ifdef NC_DOUBLE
281	ierr=NF_DEF_VAR(nid,"soildepth",NF_DOUBLE,1,idim_nsoilmx,nvarid)
282	#else
283	ierr=NF_DEF_VAR(nid,"soildepth",NF_FLOAT,1,idim_nsoilmx,nvarid)
284	#endif
285	ierr=NF_PUT_ATT_TEXT (nid,nvarid,"long_name", 20,
286	. "Soil mid-layer depth")
287	ierr=NF_PUT_ATT_TEXT (nid,nvarid,"units",1,"m")
288	ierr=NF_PUT_ATT_TEXT (nid,nvarid,"positive",4,"down")
289	ierr=NF_ENDDEF(nid) ! Leave NetCDF define mode
290	! write variable
291	#ifdef NC_DOUBLE
292	ierr=NF_PUT_VAR_DOUBLE (nid,nvarid,mlayer)
293	#else
294	ierr=NF_PUT_VAR_REAL (nid,nvarid,mlayer)
295	#endif
296
297	c
298	c --------------------------
299	c Mesh area and conversion coefficients cov. <-> contra. <--> natural
300
301	! id(1)=idim_rlonu
302	! id(2)=idim_rlatu
303	c
304	! ierr = NF_REDEF (nid)
305	!#ifdef NC_DOUBLE
306	! ierr = NF_DEF_VAR (nid, "cu", NF_DOUBLE, 2, id,nvarid)
307	!#else
308	! ierr = NF_DEF_VAR (nid, "cu", NF_FLOAT, 2, id,nvarid)
309	!#endif
310	! ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 40,
311	! . "Conversion coefficients cov <--> natural")
312	! ierr = NF_ENDDEF(nid)
313	!#ifdef NC_DOUBLE
314	! ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,cu)
315	!#else
316	! ierr = NF_PUT_VAR_REAL (nid,nvarid,cu)
317	!#endif
318	c
319	! id(1)=idim_rlonv
320	! id(2)=idim_rlatv
321	c
322	c --------------------------
323	! ierr = NF_REDEF (nid)
324	!#ifdef NC_DOUBLE
325	! ierr = NF_DEF_VAR (nid, "cv", NF_DOUBLE, 2, id,nvarid)
326	!#else
327	! ierr = NF_DEF_VAR (nid, "cv", NF_FLOAT, 2, id,nvarid)
328	!#endif
329	! ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 40,
330	! . "Conversion coefficients cov <--> natural")
331	! ierr = NF_ENDDEF(nid)
332	!#ifdef NC_DOUBLE
333	! ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,cv)
334	!#else
335	! ierr = NF_PUT_VAR_REAL (nid,nvarid,cv)
336	!#endif
337	c
338	id(1)=idim_rlonv
339	id(2)=idim_rlatu
340	c
341	c --------------------------
342	ierr = NF_REDEF (nid)
343	#ifdef NC_DOUBLE
344	ierr = NF_DEF_VAR (nid, "aire", NF_DOUBLE, 2, id,nvarid)
345	#else
346	ierr = NF_DEF_VAR (nid, "aire", NF_FLOAT, 2, id,nvarid)
347	#endif
348	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 9,
349	. "Mesh area")
350	ierr = NF_ENDDEF(nid)
351	#ifdef NC_DOUBLE
352	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,area)
353	#else
354	ierr = NF_PUT_VAR_REAL (nid,nvarid,area)
355	#endif
356	c
357	c Surface geopotential
358	id(1)=idim_rlonv
359	id(2)=idim_rlatu
360	c
361	ierr = NF_REDEF (nid)
362	#ifdef NC_DOUBLE
363	ierr = NF_DEF_VAR (nid, "phisinit", NF_DOUBLE, 2, id,nvarid)
364	#else
365	ierr = NF_DEF_VAR (nid, "phisinit", NF_FLOAT, 2, id,nvarid)
366	#endif
367	ierr = NF_PUT_ATT_TEXT (nid, nvarid, "title", 27,
368	. "Geopotential at the surface")
369	ierr = NF_ENDDEF(nid)
370	#ifdef NC_DOUBLE
371	ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,phis)
372	#else
373	ierr = NF_PUT_VAR_REAL (nid,nvarid,phis)
374	#endif
375	c
376
377	write(,)'iniwrite: nbp_lon,nbp_lat,nbp_lev,idayref',
378	& nbp_lon,nbp_lat,nbp_lev,idayref
379	write(,)'iniwrite: rad,omeg,g,mugaz,rcp',
380	& rad,omeg,g,mugaz,rcp
381	write(,)'iniwrite: daysec,dtphys',daysec,dtphys
382
383	END

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

Download in other formats:

Original Format