Context Navigation

readaerosol.F90 @ 1237

Last change on this file since 1237 was 1223, checked in by jghattas, 15 years ago

Tout les especes aerosols pour un decenie moyenne sont maintenant regroupe dans un fichier aerosolsXXXX.nc. XXXX signifie le decenie en question.

Ce modif est fait pour ne pas avoir besoin de manipuler trop des fichiers quand tout les scenarios seront fait.

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

File size: 22.4 KB

Rev	Line
[1179]	1	! $Id: readaerosol.F90 1223 2009-08-07 14:59:33Z fairhead $
[524]	2	!
[1179]	3	MODULE readaerosol_mod
	4
	5	REAL, SAVE :: not_valid=-333.
	6
	7	CONTAINS
	8
	9	SUBROUTINE readaerosol(name_aero, type, iyr_in, klev_src, pt_ap, pt_b, pt_out, psurf, load)
	10
	11	!****************************************************************************************
	12	! This routine will read the aersosol from file.
[524]	13	!
[1179]	14	! Read a year data with get_aero_fromfile depending on aer_type :
	15	! - actuel : read year 1980
	16	! - preind : read natural data
	17	! - scenario : read one or two years and do eventually linare time interpolation
[1150]	18	!
[1179]	19	! Return pointer, pt_out, to the year read or result from interpolation
	20	!****************************************************************************************
	21	USE dimphy
[524]	22
[1179]	23	IMPLICIT NONE
[766]	24
[1179]	25	INCLUDE "iniprint.h"
[1143]	26
[1179]	27	! Input arguments
	28	CHARACTER(len=7), INTENT(IN) :: name_aero
	29	CHARACTER(len=*), INTENT(IN) :: type ! correspond to aer_type in clesphys.h
	30	INTEGER, INTENT(IN) :: iyr_in
[1150]	31
[1179]	32	! Output
	33	INTEGER, INTENT(OUT) :: klev_src
	34	REAL, POINTER, DIMENSION(:) :: pt_ap ! Pointer for describing the vertical levels
	35	REAL, POINTER, DIMENSION(:) :: pt_b ! Pointer for describing the vertical levels
	36	REAL, POINTER, DIMENSION(:,:,:) :: pt_out ! The massvar distributions, DIMENSION(klon, klev_src, 12)
	37	REAL, DIMENSION(klon,12), INTENT(OUT) :: psurf ! Surface pression for 12 months
	38	REAL, DIMENSION(klon,12), INTENT(OUT) :: load ! Aerosol mass load in each column for 12 months
[1150]	39
[1179]	40	! Local variables
	41	CHARACTER(len=4) :: cyear
	42	REAL, POINTER, DIMENSION(:,:,:) :: pt_2
	43	REAL, DIMENSION(klon,12) :: psurf2, load2
	44	REAL :: p0 ! Reference pressure
	45	INTEGER :: iyr1, iyr2, klev_src2
	46	INTEGER :: it, k, i
	47	LOGICAL, PARAMETER :: lonlyone=.FALSE.
[1150]	48
[1179]	49	!****************************************************************************************
	50	! Read data depending on aer_type
	51	!
	52	!****************************************************************************************
[524]	53
[1179]	54	IF (type == 'actuel') THEN
	55	! Read and return data for year 1980
	56	!****************************************************************************************
	57	cyear='1980'
	58	! get_aero_fromfile returns pt_out allocated and initialized with data for 12 month
	59	! pt_out has dimensions (klon, klev_src, 12)
	60	CALL get_aero_fromfile(name_aero, cyear, klev_src, pt_ap, pt_b, p0, pt_out, psurf, load)
	61
[640]	62
[1179]	63	ELSE IF (type == 'preind') THEN
	64	! Read and return data from file with suffix .nat
	65	!****************************************************************************************
	66	cyear='.nat'
	67	! get_aero_fromfile returns pt_out allocated and initialized with data for 12 month
	68	! pt_out has dimensions (klon, klev_src, 12)
	69	CALL get_aero_fromfile(name_aero, cyear, klev_src, pt_ap, pt_b, p0, pt_out, psurf, load)
	70
	71	ELSE IF (type == 'scenario') THEN
	72	! Read data depending on actual year and interpolate if necessary
	73	!****************************************************************************************
	74	IF (iyr_in .LT. 1850) THEN
	75	cyear='.nat'
	76	WRITE(lunout,*) 'get_aero 1 iyr_in=', iyr_in,' ',cyear
	77	! get_aero_fromfile returns pt_out allocated and initialized with data for 12 month
	78	! pt_out has dimensions (klon, klev_src, 12)
	79	CALL get_aero_fromfile(name_aero, cyear, klev_src, pt_ap, pt_b, p0, pt_out, psurf, load)
	80
	81	ELSE IF (iyr_in .GE. 2100) THEN
	82	cyear='2100'
	83	WRITE(lunout,*) 'get_aero 2 iyr_in=', iyr_in,' ',cyear
	84	! get_aero_fromfile returns pt_out allocated and initialized with data for 12 month
	85	! pt_out has dimensions (klon, klev_src, 12)
	86	CALL get_aero_fromfile(name_aero, cyear, klev_src, pt_ap, pt_b, p0, pt_out, psurf, load)
	87
	88	ELSE
	89	! Read data from 2 decades and interpolate to actual year
	90	! a) from actual 10-yr-period
	91	IF (iyr_in.LT.1900) THEN
	92	iyr1 = 1850
	93	iyr2 = 1900
	94	ELSE IF (iyr_in.GE.1900.AND.iyr_in.LT.1920) THEN
	95	iyr1 = 1900
	96	iyr2 = 1920
	97	ELSE
	98	iyr1 = INT(iyr_in/10)*10
	99	iyr2 = INT(1+iyr_in/10)*10
	100	ENDIF
	101
	102	WRITE(cyear,'(I4)') iyr1
	103	WRITE(lunout,*) 'get_aero 3 iyr_in=', iyr_in,' ',cyear
	104	! get_aero_fromfile returns pt_out allocated and initialized with data for 12 month
	105	! pt_out has dimensions (klon, klev_src, 12)
	106	CALL get_aero_fromfile(name_aero, cyear, klev_src, pt_ap, pt_b, p0, pt_out, psurf, load)
	107
	108	! If to read two decades:
	109	IF (.NOT.lonlyone) THEN
	110
	111	! b) from the next following one
	112	WRITE(cyear,'(I4)') iyr2
	113	WRITE(lunout,*) 'get_aero 4 iyr_in=', iyr_in,' ',cyear
	114
	115	NULLIFY(pt_2)
	116	! get_aero_fromfile returns pt_2 allocated and initialized with data for 12 month
	117	! pt_2 has dimensions (klon, klev_src, 12)
	118	CALL get_aero_fromfile(name_aero, cyear, klev_src2, pt_ap, pt_b, p0, pt_2, psurf2, load2)
	119	! Test for same number of vertical levels
	120	IF (klev_src /= klev_src2) THEN
	121	WRITE(lunout,*) 'Two aerosols files with different number of vertical levels is not allowded'
	122	CALL abort_gcm('readaersosol','Error in number of vertical levels',1)
	123	END IF
	124
	125	! Linare interpolate to the actual year:
	126	DO it=1,12
	127	DO k=1,klev_src
	128	DO i = 1, klon
	129	pt_out(i,k,it) = &
	130	pt_out(i,k,it) - FLOAT(iyr_in-iyr1)/FLOAT(iyr2-iyr1) * &
	131	(pt_out(i,k,it) - pt_2(i,k,it))
	132	END DO
	133	END DO
[1143]	134
[1179]	135	DO i = 1, klon
	136	psurf(i,it) = &
	137	psurf(i,it) - FLOAT(iyr_in-iyr1)/FLOAT(iyr2-iyr1) * &
	138	(psurf(i,it) - psurf2(i,it))
[640]	139
[1179]	140	load(i,it) = &
	141	load(i,it) - FLOAT(iyr_in-iyr1)/FLOAT(iyr2-iyr1) * &
	142	(load(i,it) - load2(i,it))
	143	END DO
	144	END DO
[524]	145
[1179]	146	! Deallocate pt_2 no more needed
	147	DEALLOCATE(pt_2)
	148
	149	END IF ! lonlyone
	150	END IF ! iyr_in .LT. 1850
[524]	151
[1179]	152	ELSE
	153	WRITE(lunout,*)'This option is not implemented : aer_type = ', type
	154	CALL abort_gcm('readaerosol','Error : aer_type parameter not accepted',1)
	155	END IF ! type
[640]	156
[1151]	157
[1179]	158	END SUBROUTINE readaerosol
[1143]	159
	160
[1179]	161	SUBROUTINE get_aero_fromfile(varname, cyr, klev_src, pt_ap, pt_b, p0, pt_year, psurf_out, load_out)
	162	!****************************************************************************************
	163	! Read 12 month aerosol from file and distribute to local process on physical grid.
	164	! Vertical levels, klev_src, may differ from model levels if new file format.
	165	!
	166	! For mpi_root and master thread :
	167	! 1) Open file
	168	! 2) Find vertical dimension klev_src
	169	! 3) Read field month by month
	170	! 4) Close file
	171	! 5) Transform the global field from 2D(iim, jjp+1) to 1D(klon_glo)
	172	! - Also the levels and the latitudes have to be inversed
	173	!
	174	! For all processes and threads :
	175	! 6) Scatter global field(klon_glo) to local process domain(klon)
	176	! 7) Test for negative values
	177	!****************************************************************************************
[1150]	178
[1179]	179	USE netcdf
	180	USE dimphy
	181	USE mod_grid_phy_lmdz
	182	USE mod_phys_lmdz_para
[1183]	183	USE iophy, ONLY : io_lon, io_lat
[1179]	184
	185	IMPLICIT NONE
[524]	186
[1179]	187	INCLUDE "dimensions.h"
	188	INCLUDE "iniprint.h"
[524]	189
[1179]	190	! Input argumets
	191	CHARACTER(len=7), INTENT(IN) :: varname
	192	CHARACTER(len=4), INTENT(IN) :: cyr
[1150]	193
[1179]	194	! Output arguments
	195	INTEGER, INTENT(OUT) :: klev_src ! Number of vertical levels in file
	196	REAL, POINTER, DIMENSION(:) :: pt_ap ! Pointer for describing the vertical levels
	197	REAL, POINTER, DIMENSION(:) :: pt_b ! Pointer for describing the vertical levels
	198	REAL :: p0 ! Reference pressure value
	199	REAL, POINTER, DIMENSION(:,:,:) :: pt_year ! Pointer-variabale from file, 12 month, grid : klon,klev_src
	200	REAL, DIMENSION(klon,12), INTENT(OUT) :: psurf_out ! Surface pression for 12 months
	201	REAL, DIMENSION(klon,12), INTENT(OUT) :: load_out ! Aerosol mass load in each column
[1150]	202
[1179]	203	! Local variables
	204	CHARACTER(len=30) :: fname
[1223]	205	CHARACTER(len=8) :: filename='aerosols'
[1179]	206	CHARACTER(len=30) :: cvar
	207	INTEGER :: ncid, dimid, varid
	208	INTEGER :: imth, i, j, k, ierr
	209	REAL :: npole, spole
	210	REAL, ALLOCATABLE, DIMENSION(:,:,:) :: varmth
	211	REAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: varyear ! Global variable read from file, 12 month
	212	REAL, ALLOCATABLE, DIMENSION(:,:,:) :: varyear_glo1D !(klon_glo, klev_src, 12)
	213	REAL, ALLOCATABLE, DIMENSION(:) :: varktmp
[524]	214
[1179]	215	REAL, DIMENSION(iim,jjm+1,12) :: psurf_glo2D ! Surface pression for 12 months on dynamics global grid
	216	REAL, DIMENSION(klon_glo,12) :: psurf_glo1D ! -"- on physical global grid
	217	REAL, DIMENSION(iim,jjm+1,12) :: load_glo2D ! Load for 12 months on dynamics global grid
	218	REAL, DIMENSION(klon_glo,12) :: load_glo1D ! -"- on physical global grid
	219	REAL, DIMENSION(iim,jjm+1) :: vartmp
[1183]	220	REAL, DIMENSION(iim) :: lon_src ! longitudes in file
	221	REAL, DIMENSION(jjm+1) :: lat_src, lat_src_inv ! latitudes in file
	222	LOGICAL :: new_file ! true if new file format detected
	223	LOGICAL :: invert_lat ! true if the field has to be inverted for latitudes
[766]	224
[524]	225
[1179]	226	! Deallocate pointers
	227	IF (ASSOCIATED(pt_ap)) DEALLOCATE(pt_ap)
	228	IF (ASSOCIATED(pt_b)) DEALLOCATE(pt_b)
[766]	229
[1179]	230	!$OMP MASTER
	231	IF (is_mpi_root) THEN
[524]	232
[1179]	233	! 1) Open file
	234	!****************************************************************************************
[1223]	235	fname = filename//cyr//'.nc'
[1150]	236
[1179]	237	WRITE(lunout,*) 'reading ', TRIM(fname)
	238	CALL check_err( nf90_open(TRIM(fname), NF90_NOWRITE, ncid) )
[1150]	239
[1183]	240	! Test for equal longitudes and latitudes in file and model
	241	!****************************************************************************************
	242	! Read and test longitudes
	243	CALL check_err( nf90_inq_varid(ncid, 'lon', varid) )
	244	CALL check_err( nf90_get_var(ncid, varid, lon_src(:)) )
	245
[1202]	246	IF (maxval(ABS(lon_src - io_lon)) > 0.001) THEN
[1183]	247	WRITE(lunout,*) 'Problem in longitudes read from file : ',TRIM(fname)
	248	WRITE(lunout,*) 'longitudes in file ', TRIM(fname),' : ', lon_src
	249	WRITE(lunout,*) 'longitudes in model :', io_lon
	250
	251	CALL abort_gcm('get_aero_fromfile', 'longitudes are not the same in file and model',1)
	252	END IF
	253
	254	! Read and test latitudes
	255	CALL check_err( nf90_inq_varid(ncid, 'lat', varid) )
	256	CALL check_err( nf90_get_var(ncid, varid, lat_src(:)) )
	257
	258	! Invert source latitudes
	259	DO j = 1, jjm+1
	260	lat_src_inv(j) = lat_src(jjm+1 +1 -j)
	261	END DO
	262
[1202]	263	IF (maxval(ABS(lat_src - io_lat)) < 0.001) THEN
[1183]	264	! Latitudes are the same
	265	invert_lat=.FALSE.
[1202]	266	ELSE IF (maxval(ABS(lat_src_inv - io_lat)) < 0.001) THEN
[1183]	267	! Inverted source latitudes correspond to model latitudes
	268	WRITE(lunout,*) 'latitudes will be inverted for file : ',TRIM(fname)
	269	invert_lat=.TRUE.
	270	ELSE
	271	WRITE(lunout,*) 'Problem in latitudes read from file : ',TRIM(fname)
	272	WRITE(lunout,*) 'latitudes in file ', TRIM(fname),' : ', lat_src
	273	WRITE(lunout,*) 'latitudes in model :', io_lat
	274	CALL abort_gcm('get_aero_fromfile', 'latitudes do not correspond between file and model',1)
	275	END IF
	276
[1179]	277	! 2) Check if old or new file is avalabale.
	278	! New type of file should contain the dimension 'lev'
	279	! Old type of file should contain the dimension 'PRESNIVS'
	280	!****************************************************************************************
	281	ierr = nf90_inq_dimid(ncid, 'lev', dimid)
	282	IF (ierr /= NF90_NOERR) THEN
	283	! Coordinate axe lev not found. Check for presnivs.
	284	ierr = nf90_inq_dimid(ncid, 'PRESNIVS', dimid)
	285	IF (ierr /= NF90_NOERR) THEN
	286	! Dimension PRESNIVS not found either
	287	CALL abort_gcm('get_aero_fromfile', 'dimension lev or presnivs not in file',1)
	288	ELSE
	289	! Old file found
	290	new_file=.FALSE.
	291	WRITE(lunout,*) 'Vertical interpolation for ',TRIM(varname),' will not be done'
	292	END IF
	293	ELSE
	294	! New file found
	295	new_file=.TRUE.
	296	WRITE(lunout,*) 'Vertical interpolation for ',TRIM(varname),' will be done'
	297	END IF
	298
	299	! 2) Find vertical dimension klev_src
	300	!****************************************************************************************
	301	CALL check_err( nf90_inquire_dimension(ncid, dimid, len = klev_src) )
	302
	303	! Allocate variables depending on the number of vertical levels
	304	ALLOCATE(varmth(iim, jjm+1, klev_src), varyear(iim, jjm+1, klev_src, 12), stat=ierr)
	305	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 1',1)
[1150]	306
[1179]	307	ALLOCATE(pt_ap(klev_src), pt_b(klev_src), varktmp(klev_src), stat=ierr)
	308	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 2',1)
[1150]	309
[1179]	310	! 3) Read all variables from file
	311	! There is 2 options for the file structure :
	312	! new_file=TRUE : read varyear, ps, pt_ap and pt_b
	313	! new_file=FALSE : read varyear month by month
	314	!****************************************************************************************
[1150]	315
[1179]	316	IF (new_file) THEN
[1150]	317
[1179]	318	! ++) Read the aerosol concentration month by month and concatenate to total variable varyear
	319	!****************************************************************************************
	320	! Get variable id
	321	CALL check_err( nf90_inq_varid(ncid, TRIM(varname), varid) )
	322
	323	! Get the variable
	324	CALL check_err( nf90_get_var(ncid, varid, varyear(:,:,:,:)) )
	325
	326	! ++) Read surface pression, 12 month in one variable
	327	!****************************************************************************************
	328	! Get variable id
	329	CALL check_err( nf90_inq_varid(ncid, "ps", varid) )
	330	! Get the variable
	331	CALL check_err( nf90_get_var(ncid, varid, psurf_glo2D) )
	332
	333	! ++) Read mass load, 12 month in one variable
	334	!****************************************************************************************
	335	! Get variable id
	336	CALL check_err( nf90_inq_varid(ncid, "load_"//TRIM(varname), varid) )
	337	! Get the variable
	338	CALL check_err( nf90_get_var(ncid, varid, load_glo2D) )
	339
	340	! ++) Read ap
	341	!****************************************************************************************
	342	! Get variable id
	343	CALL check_err( nf90_inq_varid(ncid, "ap", varid) )
	344	! Get the variable
	345	CALL check_err( nf90_get_var(ncid, varid, pt_ap) )
[1150]	346
[1179]	347	! ++) Read b
	348	!****************************************************************************************
	349	! Get variable id
	350	CALL check_err( nf90_inq_varid(ncid, "b", varid) )
	351	! Get the variable
	352	CALL check_err( nf90_get_var(ncid, varid, pt_b) )
[1150]	353
[1179]	354	! ++) Read p0 : reference pressure
	355	!****************************************************************************************
	356	! Get variable id
	357	CALL check_err( nf90_inq_varid(ncid, "p0", varid) )
	358	! Get the variable
	359	CALL check_err( nf90_get_var(ncid, varid, p0) )
[1150]	360
[524]	361
[1179]	362	ELSE ! old file
[524]	363
[1179]	364	! ++) Read the aerosol concentration month by month and concatenate to total variable varyear
	365	!****************************************************************************************
	366	DO imth=1, 12
	367	IF (imth.EQ.1) THEN
	368	cvar=TRIM(varname)//'JAN'
	369	ELSE IF (imth.EQ.2) THEN
	370	cvar=TRIM(varname)//'FEB'
	371	ELSE IF (imth.EQ.3) THEN
	372	cvar=TRIM(varname)//'MAR'
	373	ELSE IF (imth.EQ.4) THEN
	374	cvar=TRIM(varname)//'APR'
	375	ELSE IF (imth.EQ.5) THEN
	376	cvar=TRIM(varname)//'MAY'
	377	ELSE IF (imth.EQ.6) THEN
	378	cvar=TRIM(varname)//'JUN'
	379	ELSE IF (imth.EQ.7) THEN
	380	cvar=TRIM(varname)//'JUL'
	381	ELSE IF (imth.EQ.8) THEN
	382	cvar=TRIM(varname)//'AUG'
	383	ELSE IF (imth.EQ.9) THEN
	384	cvar=TRIM(varname)//'SEP'
	385	ELSE IF (imth.EQ.10) THEN
	386	cvar=TRIM(varname)//'OCT'
	387	ELSE IF (imth.EQ.11) THEN
	388	cvar=TRIM(varname)//'NOV'
	389	ELSE IF (imth.EQ.12) THEN
	390	cvar=TRIM(varname)//'DEC'
	391	END IF
	392
	393	! Get variable id
	394	CALL check_err( nf90_inq_varid(ncid, TRIM(cvar), varid) )
	395
	396	! Get the variable
	397	CALL check_err( nf90_get_var(ncid, varid, varmth) )
	398
	399	! Store in variable for the whole year
	400	varyear(:,:,:,imth)=varmth(:,:,:)
	401
	402	END DO
[1150]	403
[1179]	404	! Putting dummy
	405	psurf_glo2D(:,:,:) = not_valid
	406	load_glo2D(:,:,:) = not_valid
	407	pt_ap(:) = not_valid
	408	pt_b(:) = not_valid
[524]	409
[1179]	410	END IF
[524]	411
[1179]	412	! 4) Close file
	413	!****************************************************************************************
	414	CALL check_err( nf90_close(ncid) )
	415
[524]	416
[1179]	417	! 5) Transform the global field from 2D(iim, jjp+1) to 1D(klon_glo)
	418	!****************************************************************************************
	419	! Test if vertical levels have to be inversed
[782]	420
[1179]	421	IF ((pt_b(1) < pt_b(klev_src)) .OR. .NOT. new_file) THEN
[1183]	422	WRITE(lunout,*) 'Vertical axis in file ',TRIM(fname), ' needs to be inverted'
[1179]	423	WRITE(lunout,*) 'before pt_ap = ', pt_ap
	424	WRITE(lunout,*) 'before pt_b = ', pt_b
	425
	426	! Inverse vertical levels for varyear
	427	DO imth=1, 12
	428	varmth(:,:,:) = varyear(:,:,:,imth) ! use varmth temporarly
	429	DO k=1, klev_src
	430	DO j=1, jjm+1
	431	DO i=1,iim
	432	varyear(i,j,k,imth) = varmth(i,j,klev_src+1-k)
	433	END DO
	434	END DO
	435	END DO
	436	END DO
	437
	438	! Inverte vertical axes for pt_ap and pt_b
	439	varktmp(:) = pt_ap(:)
	440	DO k=1, klev_src
	441	pt_ap(k) = varktmp(klev_src+1-k)
	442	END DO
[782]	443
[1179]	444	varktmp(:) = pt_b(:)
	445	DO k=1, klev_src
	446	pt_b(k) = varktmp(klev_src+1-k)
	447	END DO
	448	WRITE(lunout,*) 'after pt_ap = ', pt_ap
	449	WRITE(lunout,*) 'after pt_b = ', pt_b
[524]	450
[1179]	451	ELSE
	452	WRITE(lunout,*) 'Vertical axis in file ',TRIM(fname), ' is ok, no vertical inversion is done'
	453	WRITE(lunout,*) 'pt_ap = ', pt_ap
	454	WRITE(lunout,*) 'pt_b = ', pt_b
	455	END IF
[524]	456
[1183]	457	! - Invert latitudes if necessary
[1179]	458	DO imth=1, 12
[1183]	459	IF (invert_lat) THEN
	460
	461	! Invert latitudes for the variable
	462	varmth(:,:,:) = varyear(:,:,:,imth) ! use varmth temporarly
	463	DO k=1,klev_src
	464	DO j=1,jjm+1
	465	DO i=1,iim
	466	varyear(i,j,k,imth) = varmth(i,jjm+1+1-j,k)
	467	END DO
	468	END DO
	469	END DO
	470
	471	! Invert latitudes for surface pressure
	472	vartmp(:,:) = psurf_glo2D(:,:,imth)
	473	DO j=1, jjm+1
[1179]	474	DO i=1,iim
[1183]	475	psurf_glo2D(i,j,imth)= vartmp(i,jjm+1+1-j)
[1179]	476	END DO
	477	END DO
[1183]	478
	479	! Invert latitudes for the load
	480	vartmp(:,:) = load_glo2D(:,:,imth)
	481	DO j=1, jjm+1
	482	DO i=1,iim
	483	load_glo2D(i,j,imth)= vartmp(i,jjm+1+1-j)
	484	END DO
[1179]	485	END DO
[1183]	486	END IF ! invert_lat
	487
[1179]	488	! Do zonal mead at poles and distribut at whole first and last latitude
	489	DO k=1, klev_src
	490	npole=0. ! North pole, j=1
	491	spole=0. ! South pole, j=jjm+1
	492	DO i=1,iim
	493	npole = npole + varyear(i,1,k,imth)
	494	spole = spole + varyear(i,jjm+1,k,imth)
	495	END DO
	496	npole = npole/FLOAT(iim)
	497	spole = spole/FLOAT(iim)
	498	varyear(:,1, k,imth) = npole
	499	varyear(:,jjm+1,k,imth) = spole
	500	END DO
	501	END DO ! imth
[1183]	502
[1179]	503	ALLOCATE(varyear_glo1D(klon_glo, klev_src, 12), stat=ierr)
	504	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 3',1)
	505
	506	! Transform from 2D to 1D field
	507	CALL grid2Dto1D_glo(varyear,varyear_glo1D)
	508	CALL grid2Dto1D_glo(psurf_glo2D,psurf_glo1D)
	509	CALL grid2Dto1D_glo(load_glo2D,load_glo1D)
	510
	511	END IF ! is_mpi_root
[1189]	512	!$OMP END MASTER
	513	!$OMP BARRIER
[1179]	514
	515	! 6) Distribute to all processes
	516	! Scatter global field(klon_glo) to local process domain(klon)
	517	! and distribute klev_src to all processes
	518	!****************************************************************************************
[524]	519
[1179]	520	! Distribute klev_src
	521	CALL bcast(klev_src)
[524]	522
[1179]	523	! Allocate and distribute pt_ap and pt_b
	524	IF (.NOT. ASSOCIATED(pt_ap)) THEN ! if pt_ap is allocated also pt_b is allocated
	525	ALLOCATE(pt_ap(klev_src), pt_b(klev_src), stat=ierr)
	526	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 4',1)
	527	END IF
	528	CALL bcast(pt_ap)
	529	CALL bcast(pt_b)
[524]	530
[1179]	531	! Allocate space for output pointer variable at local process
	532	IF (ASSOCIATED(pt_year)) DEALLOCATE(pt_year)
	533	ALLOCATE(pt_year(klon, klev_src, 12), stat=ierr)
	534	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 5',1)
[524]	535
[1179]	536	! Scatter global field to local domain at local process
	537	CALL scatter(varyear_glo1D, pt_year)
	538	CALL scatter(psurf_glo1D, psurf_out)
	539	CALL scatter(load_glo1D, load_out)
[524]	540
[1179]	541	! 7) Test for negative values
	542	!****************************************************************************************
	543	IF (MINVAL(pt_year) < 0.) THEN
	544	WRITE(lunout,*) 'Warning! Negative values read from file :', fname
	545	END IF
[524]	546
[1179]	547	END SUBROUTINE get_aero_fromfile
[524]	548
	549
[1179]	550	SUBROUTINE check_err(status)
	551	USE netcdf
	552	IMPLICIT NONE
[524]	553
[1179]	554	INCLUDE "iniprint.h"
	555	INTEGER, INTENT (IN) :: status
[524]	556
[1179]	557	IF (status /= NF90_NOERR) THEN
	558	WRITE(lunout,*) 'Error in get_aero_fromfile ',status
	559	CALL abort_gcm('get_aero_fromfile',trim(nf90_strerror(status)),1)
	560	END IF
[524]	561
[1179]	562	END SUBROUTINE check_err
	563
	564
	565	END MODULE readaerosol_mod

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

Download in other formats:

Original Format