Context Navigation

readaerosol.F90 @ 1157

Last change on this file since 1157 was 1151, checked in by jghattas, 17 years ago

readaerosol : ENDIF mal place. Bug existant depuis la merge avec la branche LMDZ4_V3_patches. Le bug existe sur cette branche.

etat0_netcdf : manque des variables dans call conf_phys.

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 24.0 KB

Line
1	!
2	! $Header$
3	!
4	SUBROUTINE readaerosol (id_aero, r_day, first, massvar_p)
5
6	USE dimphy, ONLY : klev
7	USE mod_grid_phy_lmdz, klon=>klon_glo
8	USE mod_phys_lmdz_para
9
10	IMPLICIT NONE
11
12	! Content:
13	! --------
14	! This routine reads in monthly mean values of massvar aerosols and
15	! returns a linearly interpolated dayly-mean field.
16	!
17	!
18	! Author:
19	! -------
20	! Johannes Quaas (quaas@lmd.jussieu.fr)
21	! Celine Deandreis & Anne Cozic
22	! 19/02/09
23	!
24
25	!
26	! Include-files:
27	! --------------
28	INCLUDE "YOMCST.h"
29	INCLUDE "chem.h"
30	INCLUDE "dimensions.h"
31	INCLUDE "temps.h"
32	INCLUDE "clesphys.h"
33	INCLUDE "iniprint.h"
34	!
35	! Input:
36	! ------
37	INTEGER, INTENT(in) :: id_aero
38	REAL, INTENT(in) :: r_day ! Day of integration
39	LOGICAL, INTENT(in) :: first ! First timestep
40	! (and therefore initialization necessary)
41	!
42	! Output:
43	! -------
44	REAL, INTENT(out) :: massvar_p(klon_omp,klev) ! Mass of massvar (monthly mean data,from file) [ug AIBCM/m3]
45
46	!
47	! Local Variables:
48	! ----------------
49	INTEGER :: i, ig, k, it
50	INTEGER :: j, iday, iyr, iyr1, iyr2
51	INTEGER :: im, day1, day2, im2
52	INTEGER, PARAMETER :: ny=jjm+1
53	CHARACTER(len=4) :: cyear
54	CHARACTER(len=7),DIMENSION(8) :: name_aero
55	REAL :: var_1(iim, jjm+1, klev, 12)
56	REAL, DIMENSION(klon,klev) :: massvar
57	REAL, DIMENSION(iim, jjm+1, klev, 12) :: var_2 ! The massvar distributions
58	REAL, ALLOCATABLE, DIMENSION(:,:,:,:), SAVE :: var ! VAR in right dimension
59	REAL, ALLOCATABLE, DIMENSION(:,:,:),SAVE :: var_out
60	!$OMP THREADPRIVATE(var,var_out)
61
62	LOGICAL :: lnewday
63	LOGICAL, PARAMETER :: lonlyone=.FALSE.
64	LOGICAL,SAVE :: first2=.TRUE.
65	!$OMP THREADPRIVATE(first2)
66
67	! Variable pour definir a partir d'un numero d'aerosol, son nom
68	name_aero(1) = "SSSSM "
69	name_aero(2) = "ASSSM "
70	name_aero(3) = "ASBCM "
71	name_aero(4) = "ASPOMM "
72	name_aero(5) = "SO4 "
73	name_aero(6) = "CIDUSTM"
74	name_aero(7) = "AIBCM "
75	name_aero(8) = "AIPOMM "
76
77	!$OMP MASTER
78	IF (first2) THEN
79	ALLOCATE( var(klon, klev, 12,8) )
80	ALLOCATE( var_out(klon, klev,8))
81	first2=.FALSE.
82
83	IF (aer_type /= 'actuel ' .AND. aer_type /= 'preind ' .AND. &
84	aer_type /= 'scenario') THEN
85	WRITE(lunout,)' Warning *'
86	WRITE(lunout,*)'Option aer_type non prevu pour les aerosols = ', &
87	aer_type
88	CALL abort_gcm('readaerosol','Error : aer_type parameter not accepted',1)
89	ENDIF
90	ENDIF
91
92	IF (is_mpi_root) THEN
93
94	iday = INT(r_day)
95	! Get the year of the run
96	iyr = iday/360
97	! Get the day of the actual year:
98	iday = iday-iyr*360
99	! 0.02 is about 0.5/24, namly less than half an hour
100	lnewday = (r_day-FLOAT(iday).LT.0.02)
101
102	! ---------------------------------------------
103	! All has to be done only, if a new day begins!
104	! ---------------------------------------------
105
106	IF (lnewday.OR.first) THEN
107
108	im = iday/30 +1 ! the actual month
109	! annee_ref is the initial year (defined in temps.h)
110	iyr = iyr + annee_ref
111
112	! Do I have to read new data? (Is this the first day of a year?)
113	IF (first.OR.iday.EQ.1.) THEN
114	! Initialize values
115	DO it=1,12
116	DO k=1,klev
117	DO i=1,klon
118	var(i,k,it,id_aero)=0.
119	ENDDO
120	ENDDO
121	ENDDO
122
123
124	IF (aer_type == 'actuel ') then
125
126	cyear='1980'
127	CALL get_aero_fromfile(cyear, var_1, name_aero(id_aero))
128
129	ELSE IF (aer_type == 'preind ') THEN
130
131	cyear='.nat'
132	CALL get_aero_fromfile(cyear, var_1, name_aero(id_aero))
133
134	ELSE ! aer_type=scenario
135
136	IF (iyr .LT. 1850) THEN
137	cyear='.nat'
138	WRITE(lunout,*) 'get_aero iyr=', iyr,' ',cyear
139	CALL get_aero_fromfile(cyear, var_1, name_aero(id_aero))
140
141	ELSE IF (iyr .GE. 2100) THEN
142	cyear='2100'
143	WRITE(lunout,*) 'get_aero iyr=', iyr,' ',cyear
144	CALL get_aero_fromfile(cyear, var_1, name_aero(id_aero))
145
146	ELSE ! Read data from 2 decades
147	! a) from actual 10-yr-period
148	IF (iyr.LT.1900) THEN
149	iyr1 = 1850
150	iyr2 = 1900
151	ELSE IF (iyr.GE.1900.AND.iyr.LT.1920) THEN
152	iyr1 = 1900
153	iyr2 = 1920
154	ELSE
155	iyr1 = INT(iyr/10)*10
156	iyr2 = INT(1+iyr/10)*10
157	ENDIF
158
159	WRITE(cyear,'(I4)') iyr1
160	WRITE(lunout,*) 'get_aero iyr=', iyr,' ',cyear
161	CALL get_aero_fromfile(cyear, var_1, name_aero(id_aero))
162
163	! If to read two decades:
164	IF (.NOT.lonlyone) THEN
165
166	! b) from the next following one
167	WRITE(cyear,'(I4)') iyr2
168	WRITE(lunout,*) 'get_aero iyr=', iyr,' ',cyear
169
170	CALL get_aero_fromfile(cyear, var_2, name_aero(id_aero))
171
172	! Interpolate linarily to the actual year:
173	DO it=1,12
174	DO k=1,klev
175	DO j=1,jjm
176	DO i=1,iim
177	var_1(i,j,k,it) = &
178	var_1(i,j,k,it) - FLOAT(iyr-iyr1)/FLOAT(iyr2-iyr1) * &
179	(var_1(i,j,k,it) - var_2(i,j,k,it))
180	ENDDO
181	ENDDO
182	ENDDO
183	ENDDO
184
185	ENDIF ! lonlyone
186	ENDIF ! iyr .LT. 1850
187	ENDIF ! aer_type
188
189	! Transform the horizontal 2D-field into the physics-field
190	! (Also the levels and the latitudes have to be inversed)
191
192	DO it=1,12
193	DO k=1, klev
194	! a) at the poles, use the zonal mean:
195	DO i=1,iim
196	! North pole
197	var(1,k,it,id_aero) = &
198	var(1,k,it,id_aero)+ var_1(i,jjm+1,klev+1-k,it)
199	! South pole
200	var(klon,k,it,id_aero)= &
201	var(klon,k,it,id_aero)+ var_1(i,1,klev+1-k,it)
202	ENDDO
203	var(1,k,it,id_aero) = var(1,k,it,id_aero)/FLOAT(iim)
204	var(klon,k,it,id_aero)= var(klon,k,it,id_aero)/FLOAT(iim)
205
206	! b) the values between the poles:
207	ig=1
208	DO j=2,jjm
209	DO i=1,iim
210	ig=ig+1
211
212	IF (ig.GT.klon) THEN
213	WRITE(lunout,*) 'Attention dans readaerosol', &
214	name_aero(id_aero), 'ig > klon', ig, klon
215	ENDIF
216
217	var(ig,k,it,id_aero) = var_1(i,jjm+1+1-j,klev+1-k,it)
218
219	ENDDO
220	ENDDO
221	IF (ig.NE.klon-1) THEN
222	PRINT *,'Error in readaerosol', name_aero(id_aero), 'conversion'
223	CALL abort_gcm('readaerosol','Error in readaerosol 1',1)
224	ENDIF
225	ENDDO ! Loop over k (vertical)
226	ENDDO ! Loop over it (months)
227
228	ENDIF ! Had to read new data?
229
230
231	! Interpolate to actual day:
232	IF (iday.LT.im*30-15) THEN
233	! in the first half of the month use month before and actual month
234	im2=im-1
235	day2 = im2*30-15
236	day1 = im2*30+15
237	IF (im2.LE.0) THEN
238	! the month is january, thus the month before december
239	im2=12
240	ENDIF
241	DO k=1,klev
242	DO i=1,klon
243	massvar(i,k) = &
244	var(i,k,im2,id_aero)- FLOAT(iday-day2)/FLOAT(day1-day2) * &
245	(var(i,k,im2,id_aero) - var(i,k,im,id_aero))
246
247	IF (massvar(i,k).LT.0.) THEN
248	IF (iday-day2.LT.0.) WRITE(lunout,*) 'iday-day2',iday-day2
249	IF (var(i,k,im2,id_aero) - var(i,k,im,id_aero).LT.0.) THEN
250	PRINT*, name_aero(id_aero),'(i,k,im2)- ', &
251	name_aero(id_aero),'(i,k,im)', &
252	var(i,k,im2,id_aero) - var(i,k,im,id_aero)
253	ENDIF
254
255	IF (day1-day2.LT.0.) WRITE(lunout,*) 'day1-day2',day1-day2
256	WRITE(lunout,*) 'stop',name_aero(id_aero)
257	CALL abort_gcm('readaerosol','Error in interpolation 2',1)
258	ENDIF
259	ENDDO
260	ENDDO
261	ELSE
262	! the second half of the month
263	im2=im+1
264	IF (im2.GT.12) THEN
265	! the month is december, the following thus january
266	im2=1
267	ENDIF
268	day2 = im*30-15
269	day1 = im*30+15
270	DO k=1,klev
271	DO i=1,klon
272	massvar(i,k) = &
273	var(i,k,im2,id_aero) - FLOAT(iday-day2)/FLOAT(day1-day2) * &
274	(var(i,k,im2,id_aero) - var(i,k,im,id_aero))
275
276	IF (massvar(i,k).LT.0.) THEN
277	IF (iday-day2.LT.0.) &
278	WRITE(lunout,*) 'iday-day2',iday-day2
279	IF (var(i,k,im2,id_aero) - var(i,k,im,id_aero).LT.0.) THEN
280	WRITE(lunout,*) name_aero(id_aero),'(i,k,im2)-', &
281	name_aero(id_aero),'(i,k,im)', &
282	var(i,k,im2,id_aero) - var(i,k,im,id_aero)
283	ENDIF
284	IF (day1-day2.LT.0.) &
285	WRITE(lunout,*) 'day1-day2',day1-day2
286	WRITE(lunout,*) 'stop',name_aero(id_aero)
287	CALL abort_gcm('readaerosol','Error in interpolation 3',1)
288	ENDIF
289	ENDDO
290	ENDDO
291	ENDIF
292
293
294	! The massvar concentration [molec cm-3] is read in.
295	! Convert it into mass [ug VAR/m3]
296	! masse_var in [g/mol], n_avogadro in [molec/mol]
297	! The sulfate mass [ug VAR/m3] is read in.
298	DO k=1,klev
299	DO i=1,klon
300	var_out(i,k,id_aero) = massvar(i,k)
301	IF (var_out(i,k,id_aero).LT.0) THEN
302	PRINT*, 'Attention il y a un probleme dans readaerosol', &
303	name_aero(id_aero), 'la masse est negative'
304	CALL abort_gcm('readaerosol','Error in readaerosol 4',1)
305	ENDIF
306	ENDDO
307	ENDDO
308	ELSE ! if no new day, use old data:
309	DO k=1,klev
310	DO i=1,klon
311	massvar(i,k) = var_out(i,k,id_aero)
312	IF (var_out(i,k,id_aero).LT.0) THEN
313	PRINT*, 'Attention il y a un probleme dans readaerosol', &
314	name_aero(id_aero), 'la masse est negative'
315	CALL abort_gcm('readaerosol','Error in readaerosol 5',1)
316	ENDIF
317	ENDDO
318	ENDDO
319
320
321	ENDIF ! Did I have to do anything (was it a new day?)
322
323	ENDIF ! phy_rank==0
324
325	!$OMP END MASTER
326	CALL Scatter(massvar,massvar_p)
327
328	END SUBROUTINE readaerosol
329
330
331
332
333
334	!-----------------------------------------------------------------------------
335	! Read in /calculate pre-industrial values of sulfate
336	!-----------------------------------------------------------------------------
337
338	SUBROUTINE readaerosol_preind (id_aero, r_day, first, pi_massvar_p)
339
340	USE dimphy, ONLY : klev
341	USE mod_grid_phy_lmdz, klon=>klon_glo
342	USE mod_phys_lmdz_para
343	IMPLICIT NONE
344
345	! Content:
346	! --------
347	! This routine reads in monthly mean values of massvar aerosols and
348	! returns a linearly interpolated dayly-mean field.
349	!
350	! It does so for the preindustriel values of the massvar, to a large part
351	! analogous to the routine readmassvar above.
352	!
353	! Only Pb: Variables must be saved and don t have to be overwritten!
354	!
355	! Author:
356	! -------
357	! Celine Deandreis & Anne Cozic (LSCE) 2009
358	! Johannes Quaas (quaas@lmd.jussieu.fr)
359	! 26/06/01
360	!
361	! Modifications:
362	! --------------
363	! see above
364
365
366	! Include-files:
367	! --------------
368
369	INCLUDE "YOMCST.h"
370	INCLUDE "chem.h"
371	INCLUDE "dimensions.h"
372	INCLUDE "temps.h"
373	INCLUDE "iniprint.h"
374
375	! Input:
376	! ------
377	INTEGER, INTENT(in) :: id_aero
378	REAL, INTENT(in) :: r_day ! Day of integration
379	LOGICAL, INTENT(in) :: first ! First timestep (and therefore initialization necessary)
380
381
382	!
383	! Output:
384	! -------
385	REAL, DIMENSION(klon_omp,klev), INTENT(out) :: pi_massvar_p ! Number conc. massvar (monthly mean data, from file)
386
387
388	!
389	! Local Variables:
390	! ----------------
391	INTEGER :: i, ig, k, it
392	INTEGER :: j, iday, iyr
393	INTEGER, PARAMETER :: ny=jjm+1
394	INTEGER :: im, day1, day2, im2
395
396	REAL, DIMENSION(iim, jjm+1, klev, 12) :: pi_var_1
397	REAL, DIMENSION(klon,klev) :: pi_massvar
398	REAL, ALLOCATABLE, DIMENSION(:,:,:,:), SAVE :: pi_var ! VAR in right dimension
399	REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: pi_var_out
400	!$OMP THREADPRIVATE(pi_var,pi_var_out)
401
402	CHARACTER(len=4) :: cyear
403	CHARACTER(len=7),DIMENSION(8) :: name_aero
404	LOGICAL :: lnewday
405	LOGICAL,SAVE :: first2=.TRUE.
406	!$OMP THREADPRIVATE(first2)
407
408
409	! Variable pour definir a partir d'un numero d'aerosol, son nom
410	name_aero(1) = "SSSSM "
411	name_aero(2) = "ASSSM "
412	name_aero(3) = "ASBCM "
413	name_aero(4) = "ASPOMM "
414	name_aero(5) = "SO4 "
415	name_aero(6) = "CIDUSTM"
416	name_aero(7) = "AIBCM "
417	name_aero(8) = "AIPOMM "
418
419
420	!$OMP MASTER
421	IF (first2) THEN
422	ALLOCATE( pi_var(klon, klev, 12,8) )
423	ALLOCATE( pi_var_out(klon, klev,8))
424	first2=.FALSE.
425	ENDIF
426
427	IF (is_mpi_root) THEN
428	iday = INT(r_day)
429	! Get the year of the run
430	iyr = iday/360
431	! Get the day of the actual year:
432	iday = iday-iyr*360
433	! 0.02 is about 0.5/24, namly less than half an hour
434	lnewday = (r_day-FLOAT(iday).LT.0.02)
435
436	! ---------------------------------------------
437	! All has to be done only, if a new day begins!
438	! ---------------------------------------------
439
440	IF (lnewday.OR.first) THEN
441
442	im = iday/30 +1 ! the actual month
443	! annee_ref is the initial year (defined in temps.h)
444	iyr = iyr + annee_ref
445
446	IF (first) THEN
447	cyear='.nat'
448	CALL get_aero_fromfile(cyear,pi_var_1, name_aero(id_aero))
449
450	! Transform the horizontal 2D-field into the physics-field
451	! (Also the levels and the latitudes have to be inversed)
452
453	! Initialize field
454	DO it=1,12
455	DO k=1,klev
456	DO i=1,klon
457	pi_var(i,k,it,id_aero)=0.
458	ENDDO
459	ENDDO
460	ENDDO
461
462	WRITE(lunout,*) 'preind: finished reading', FLOAT(iim)
463	DO it=1,12
464	DO k=1, klev
465	! a) at the poles, use the zonal mean:
466	DO i=1,iim
467	! North pole
468	pi_var(1,k,it,id_aero) = &
469	pi_var(1,k,it,id_aero) + pi_var_1(i,jjm+1,klev+1-k,it)
470	! South pole
471	pi_var(klon,k,it,id_aero) = &
472	pi_var(klon,k,it,id_aero) + pi_var_1(i,1,klev+1-k,it)
473	ENDDO
474	pi_var(1,k,it,id_aero) = pi_var(1,k,it,id_aero)/FLOAT(iim)
475	pi_var(klon,k,it,id_aero) = pi_var(klon,k,it,id_aero)/FLOAT(iim)
476
477	! b) the values between the poles:
478	ig=1
479	DO j=2,jjm
480	DO i=1,iim
481	ig=ig+1
482	IF (ig.GT.klon) THEN
483	WRITE(lunout,*) 'Attention dans readaerosol_preind', &
484	name_aero(id_aero), 'ig > klon', ig, klon
485	ENDIF
486	pi_var(ig,k,it,id_aero) = &
487	pi_var_1(i,jjm+1+1-j,klev+1-k,it)
488	ENDDO
489	ENDDO
490	IF (ig.NE.klon-1) THEN
491	WRITE(lunout,*) 'Error in readaerosol_preind (', name_aero(id_aero), ' conversion)'
492	CALL abort_gcm('readaerosol_preind','Error in readaerosol_preind 1',1)
493	ENDIF
494	ENDDO ! Loop over k (vertical)
495	ENDDO ! Loop over it (months)
496
497	ENDIF ! Had to read new data?
498
499
500	! Interpolate to actual day:
501	IF (iday.LT.im*30-15) THEN
502	! in the first half of the month use month before and actual month
503	im2=im-1
504	day1 = im2*30+15
505	day2 = im2*30-15
506	IF (im2.LE.0) THEN
507	! the month is january, thus the month before december
508	im2=12
509	ENDIF
510	DO k=1,klev
511	DO i=1,klon
512	pi_massvar(i,k) = &
513	pi_var(i,k,im2,id_aero) - FLOAT(iday-day2)/FLOAT(day1-day2) * &
514	(pi_var(i,k,im2,id_aero) - pi_var(i,k,im,id_aero))
515
516	IF (pi_massvar(i,k).LT.0.) THEN
517	IF (iday-day2.LT.0.) WRITE(lunout,*) 'iday-day2',iday-day2
518	IF (pi_var(i,k,im2,id_aero) - pi_var(i,k,im,id_aero).LT.0.) THEN
519	WRITE(lunout,*) 'pi_',name_aero(id_aero),'(i,k,im2) - pi_', &
520	name_aero(id_aero),'(i,k,im)', &
521	pi_var(i,k,im2,id_aero) -pi_var(i,k,im,id_aero)
522	ENDIF
523	IF (day1-day2.LT.0.)WRITE(lunout,*) 'day1-day2',day1-day2
524	PRINT *, 'pi_',name_aero(id_aero)
525	CALL abort_gcm('readaerosol_preind','Error in readaerosol_preind 2',1)
526	ENDIF
527	ENDDO
528	ENDDO
529	ELSE
530
531	! the second half of the month
532	im2=im+1
533	day1 = im*30+15
534	IF (im2.GT.12) THEN
535	! the month is december, the following thus january
536	im2=1
537	ENDIF
538	day2 = im*30-15
539
540	DO k=1,klev
541	DO i=1,klon
542	pi_massvar(i,k) = &
543	pi_var(i,k,im2,id_aero) - FLOAT(iday-day2)/FLOAT(day1-day2) * &
544	(pi_var(i,k,im2,id_aero) - pi_var(i,k,im,id_aero))
545	IF (pi_massvar(i,k).LT.0.) THEN
546	IF (iday-day2.LT.0.) WRITE(lunout,*) 'iday-day2',iday-day2
547	IF (pi_var(i,k,im2,id_aero) - pi_var(i,k,im,id_aero).LT.0.) THEN
548	WRITE(lunout,*) 'pi_', name_aero(id_aero),'(i,k,im2) - pi_',&
549	name_aero(id_aero), '(i,k,im)',&
550	pi_var(i,k,im2,id_aero) - pi_var(i,k,im,id_aero)
551	ENDIF
552	IF (day1-day2.LT.0.) &
553	WRITE(lunout,*) 'day1-day2',day1-day2
554	PRINT *, 'pi_',name_aero(id_aero)
555	CALL abort_gcm('readaerosol_preind','Error in readaerosol_preind 3',1)
556	ENDIF
557	ENDDO
558	ENDDO
559	ENDIF
560
561
562	! The massvar concentration [molec cm-3] is read in.
563	! Convert it into mass [ug VAR/m3]
564	! masse_var in [g/mol], n_avogadro in [molec/mol]
565	DO k=1,klev
566	DO i=1,klon
567	pi_var_out(i,k,id_aero) = pi_massvar(i,k)
568	ENDDO
569	ENDDO
570
571	ELSE ! If no new day, use old data:
572	DO k=1,klev
573	DO i=1,klon
574	pi_massvar(i,k) = pi_var_out(i,k,id_aero)
575	ENDDO
576	ENDDO
577	ENDIF ! Was this the beginning of a new day?
578	ENDIF ! is_mpi_root
579
580	!$OMP END MASTER
581	CALL Scatter(pi_massvar,pi_massvar_p)
582
583	END SUBROUTINE readaerosol_preind
584
585
586
587	!-----------------------------------------------------------------------------
588	! Routine for reading VAR data from files
589	!-----------------------------------------------------------------------------
590
591
592	SUBROUTINE get_aero_fromfile (cyr, var, varname)
593	USE dimphy
594	IMPLICIT NONE
595
596	INCLUDE "netcdf.inc"
597	INCLUDE "dimensions.h"
598	INCLUDE "iniprint.h"
599
600	CHARACTER(len=4), INTENT(in) :: cyr
601	REAL, DIMENSION(iim, jjm+1, klev, 12), INTENT(out) :: var
602	CHARACTER(len=*), INTENT(in) :: varname
603
604
605	CHARACTER(len=30) :: fname
606	CHARACTER(len=30) :: cvar
607	INTEGER, DIMENSION(3) :: START, COUNT
608	INTEGER :: STATUS, NCID, VARID
609	INTEGER :: imth, i, j, k
610	INTEGER, PARAMETER :: ny=jjm+1
611	REAL, DIMENSION(iim, ny, klev) :: varmth
612	!
613	!
614	!
615	fname = TRIM(varname)//'.run'//cyr//'.cdf'
616
617	WRITE(lunout,*) 'reading ', TRIM(fname)
618	STATUS = NF_OPEN (TRIM(fname), NF_NOWRITE, NCID)
619	IF (STATUS .NE. NF_NOERR) THEN
620	WRITE(lunout,*) 'err in open ',status
621	CALL abort_gcm('get_aero_fromfile','Error in opening file',1)
622	ENDIF
623	DO imth=1, 12
624	IF (imth.EQ.1) THEN
625	cvar=TRIM(varname)//'JAN'
626	ELSEIF (imth.EQ.2) THEN
627	cvar=TRIM(varname)//'FEB'
628	ELSEIF (imth.EQ.3) THEN
629	cvar=TRIM(varname)//'MAR'
630	ELSEIF (imth.EQ.4) THEN
631	cvar=TRIM(varname)//'APR'
632	ELSEIF (imth.EQ.5) THEN
633	cvar=TRIM(varname)//'MAY'
634	ELSEIF (imth.EQ.6) THEN
635	cvar=TRIM(varname)//'JUN'
636	ELSEIF (imth.EQ.7) THEN
637	cvar=TRIM(varname)//'JUL'
638	ELSEIF (imth.EQ.8) THEN
639	cvar=TRIM(varname)//'AUG'
640	ELSEIF (imth.EQ.9) THEN
641	cvar=TRIM(varname)//'SEP'
642	ELSEIF (imth.EQ.10) THEN
643	cvar=TRIM(varname)//'OCT'
644	ELSEIF (imth.EQ.11) THEN
645	cvar=TRIM(varname)//'NOV'
646	ELSEIF (imth.EQ.12) THEN
647	cvar=TRIM(varname)//'DEC'
648	ENDIF
649	start(1)=1
650	start(2)=1
651	start(3)=1
652	COUNT(1)=iim
653	COUNT(2)=ny
654	COUNT(3)=klev
655	STATUS = NF_INQ_VARID (NCID, TRIM(cvar), VARID)
656	WRITE(lunout,*) ncid,imth,TRIM(cvar), varid
657
658	IF (STATUS .NE. NF_NOERR) WRITE(lunout,*) 'err in read ',status
659
660	#ifdef NC_DOUBLE
661	status = NF_GET_VARA_DOUBLE(NCID, VARID, START, COUNT,varmth)
662	#else
663	status = NF_GET_VARA_REAL(NCID, VARID, START, COUNT, varmth)
664	#endif
665	IF (STATUS .NE. NF_NOERR) WRITE(lunout,*) 'err in read data',status
666
667	DO k=1,klev
668	DO j=1,jjm+1
669	DO i=1,iim
670	IF (varmth(i,j,k).LT.0.) THEN
671	WRITE(lunout,*) 'this is shit'
672	WRITE(lunout,*) varname,'(',i,j,k,') =',varmth(i,j,k)
673	ENDIF
674	var(i,j,k,imth)=varmth(i,j,k)
675	ENDDO
676	ENDDO
677	ENDDO
678	ENDDO
679
680	STATUS = NF_CLOSE(NCID)
681	IF (STATUS .NE. NF_NOERR) WRITE(lunout,*) 'err in closing file',status
682
683
684	END SUBROUTINE get_aero_fromfile
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700

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

Download in other formats:

Original Format