Context Navigation

readaerosol.F90 @ 1456

Last change on this file since 1456 was 1299, checked in by Laurent Fairhead, 15 years ago

Nettoyage general pour se rapprocher des normes et éviter des erreurs a la
compilation:

tous les FLOAT() sont remplacés par des REAL()
tous les STOP dans phylmd sont remplacés par des appels à abort_gcm
le common control défini dans le fichier control.h est remplacé par le module control_mod pour éviter des messages sur l'alignement des variables dans les déclarations
des $Header$ remplacés par des $Id$ pour svn

Quelques remplacements à faire ont pu m'échapper

General cleanup of the code to try and adhere to norms and to prevent some
compilation errors:

all FLOAT() instructions have been replaced by REAL() instructions
all STOP instructions in phylmd have been replaced by calls to abort_gcm
the common block control defined in the control.h file has been replaced by the control_mod to prevent compilation warnings on the alignement of declared variables
$Header$ replaced by $Id$ for svn

Some changes which should have been made might have escaped me

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

File size: 22.5 KB

Rev	Line
[1179]	1	! $Id: readaerosol.F90 1299 2010-01-20 14:27:21Z musat $
[524]	2	!
[1179]	3	MODULE readaerosol_mod
	4
	5	REAL, SAVE :: not_valid=-333.
	6
	7	CONTAINS
	8
[1270]	9	SUBROUTINE readaerosol(name_aero, type, iyr_in, klev_src, pt_ap, pt_b, pt_out, psurf, load)
[1179]	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
[1270]	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
[1270]	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'
[1270]	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)
[1179]	61
[640]	62
[1179]	63	ELSE IF (type == 'preind') THEN
	64	! Read and return data from file with suffix .nat
	65	!****************************************************************************************
	66	cyear='.nat'
[1270]	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)
[1179]	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
[1270]	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)
[1179]	80
	81	ELSE IF (iyr_in .GE. 2100) THEN
	82	cyear='2100'
	83	WRITE(lunout,*) 'get_aero 2 iyr_in=', iyr_in,' ',cyear
[1270]	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)
[1179]	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
[1270]	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)
[1179]	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)
[1270]	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)
[1179]	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:
[1270]	126	DO it=1,12
[1179]	127	DO k=1,klev_src
	128	DO i = 1, klon
	129	pt_out(i,k,it) = &
[1299]	130	pt_out(i,k,it) - REAL(iyr_in-iyr1)/REAL(iyr2-iyr1) * &
[1179]	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
[1270]	136	psurf(i,it) = &
[1299]	137	psurf(i,it) - REAL(iyr_in-iyr1)/REAL(iyr2-iyr1) * &
[1270]	138	(psurf(i,it) - psurf2(i,it))
[640]	139
[1270]	140	load(i,it) = &
[1299]	141	load(i,it) - REAL(iyr_in-iyr1)/REAL(iyr2-iyr1) * &
[1270]	142	(load(i,it) - load2(i,it))
[1179]	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
[1270]	161	SUBROUTINE get_aero_fromfile(varname, cyr, klev_src, pt_ap, pt_b, p0, pt_year, psurf_out, load_out)
[1179]	162	!****************************************************************************************
[1270]	163	! Read 12 month aerosol from file and distribute to local process on physical grid.
[1179]	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
[1270]	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
[1270]	211	REAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: varyear ! Global variable read from file, 12 month
	212	REAL, ALLOCATABLE, DIMENSION(:,:,:) :: varyear_glo1D !(klon_glo, klev_src, 12)
[1179]	213	REAL, ALLOCATABLE, DIMENSION(:) :: varktmp
[524]	214
[1270]	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
[1179]	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
[1249]	230	IF (is_mpi_root .AND. is_omp_root) THEN
[524]	231
[1179]	232	! 1) Open file
	233	!****************************************************************************************
[1223]	234	fname = filename//cyr//'.nc'
[1150]	235
[1179]	236	WRITE(lunout,*) 'reading ', TRIM(fname)
	237	CALL check_err( nf90_open(TRIM(fname), NF90_NOWRITE, ncid) )
[1150]	238
[1183]	239	! Test for equal longitudes and latitudes in file and model
	240	!****************************************************************************************
	241	! Read and test longitudes
	242	CALL check_err( nf90_inq_varid(ncid, 'lon', varid) )
	243	CALL check_err( nf90_get_var(ncid, varid, lon_src(:)) )
	244
[1202]	245	IF (maxval(ABS(lon_src - io_lon)) > 0.001) THEN
[1183]	246	WRITE(lunout,*) 'Problem in longitudes read from file : ',TRIM(fname)
	247	WRITE(lunout,*) 'longitudes in file ', TRIM(fname),' : ', lon_src
	248	WRITE(lunout,*) 'longitudes in model :', io_lon
	249
	250	CALL abort_gcm('get_aero_fromfile', 'longitudes are not the same in file and model',1)
	251	END IF
	252
	253	! Read and test latitudes
	254	CALL check_err( nf90_inq_varid(ncid, 'lat', varid) )
	255	CALL check_err( nf90_get_var(ncid, varid, lat_src(:)) )
	256
	257	! Invert source latitudes
	258	DO j = 1, jjm+1
	259	lat_src_inv(j) = lat_src(jjm+1 +1 -j)
	260	END DO
	261
[1202]	262	IF (maxval(ABS(lat_src - io_lat)) < 0.001) THEN
[1183]	263	! Latitudes are the same
	264	invert_lat=.FALSE.
[1202]	265	ELSE IF (maxval(ABS(lat_src_inv - io_lat)) < 0.001) THEN
[1183]	266	! Inverted source latitudes correspond to model latitudes
	267	WRITE(lunout,*) 'latitudes will be inverted for file : ',TRIM(fname)
	268	invert_lat=.TRUE.
	269	ELSE
	270	WRITE(lunout,*) 'Problem in latitudes read from file : ',TRIM(fname)
	271	WRITE(lunout,*) 'latitudes in file ', TRIM(fname),' : ', lat_src
	272	WRITE(lunout,*) 'latitudes in model :', io_lat
	273	CALL abort_gcm('get_aero_fromfile', 'latitudes do not correspond between file and model',1)
	274	END IF
	275
[1179]	276	! 2) Check if old or new file is avalabale.
	277	! New type of file should contain the dimension 'lev'
	278	! Old type of file should contain the dimension 'PRESNIVS'
	279	!****************************************************************************************
	280	ierr = nf90_inq_dimid(ncid, 'lev', dimid)
	281	IF (ierr /= NF90_NOERR) THEN
	282	! Coordinate axe lev not found. Check for presnivs.
	283	ierr = nf90_inq_dimid(ncid, 'PRESNIVS', dimid)
	284	IF (ierr /= NF90_NOERR) THEN
	285	! Dimension PRESNIVS not found either
	286	CALL abort_gcm('get_aero_fromfile', 'dimension lev or presnivs not in file',1)
	287	ELSE
	288	! Old file found
	289	new_file=.FALSE.
	290	WRITE(lunout,*) 'Vertical interpolation for ',TRIM(varname),' will not be done'
	291	END IF
	292	ELSE
	293	! New file found
	294	new_file=.TRUE.
	295	WRITE(lunout,*) 'Vertical interpolation for ',TRIM(varname),' will be done'
	296	END IF
	297
	298	! 2) Find vertical dimension klev_src
	299	!****************************************************************************************
	300	CALL check_err( nf90_inquire_dimension(ncid, dimid, len = klev_src) )
	301
	302	! Allocate variables depending on the number of vertical levels
[1270]	303	ALLOCATE(varmth(iim, jjm+1, klev_src), varyear(iim, jjm+1, klev_src, 12), stat=ierr)
[1179]	304	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 1',1)
[1150]	305
[1179]	306	ALLOCATE(pt_ap(klev_src), pt_b(klev_src), varktmp(klev_src), stat=ierr)
	307	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 2',1)
[1150]	308
[1179]	309	! 3) Read all variables from file
	310	! There is 2 options for the file structure :
	311	! new_file=TRUE : read varyear, ps, pt_ap and pt_b
	312	! new_file=FALSE : read varyear month by month
	313	!****************************************************************************************
[1150]	314
[1179]	315	IF (new_file) THEN
[1150]	316
[1179]	317	! ++) Read the aerosol concentration month by month and concatenate to total variable varyear
	318	!****************************************************************************************
	319	! Get variable id
	320	CALL check_err( nf90_inq_varid(ncid, TRIM(varname), varid) )
	321
	322	! Get the variable
	323	CALL check_err( nf90_get_var(ncid, varid, varyear(:,:,:,:)) )
	324
[1270]	325	! ++) Read surface pression, 12 month in one variable
[1179]	326	!****************************************************************************************
	327	! Get variable id
	328	CALL check_err( nf90_inq_varid(ncid, "ps", varid) )
	329	! Get the variable
	330	CALL check_err( nf90_get_var(ncid, varid, psurf_glo2D) )
	331
	332	! ++) Read mass load, 12 month in one variable
	333	!****************************************************************************************
	334	! Get variable id
	335	CALL check_err( nf90_inq_varid(ncid, "load_"//TRIM(varname), varid) )
	336	! Get the variable
	337	CALL check_err( nf90_get_var(ncid, varid, load_glo2D) )
	338
	339	! ++) Read ap
	340	!****************************************************************************************
	341	! Get variable id
	342	CALL check_err( nf90_inq_varid(ncid, "ap", varid) )
	343	! Get the variable
	344	CALL check_err( nf90_get_var(ncid, varid, pt_ap) )
[1150]	345
[1179]	346	! ++) Read b
	347	!****************************************************************************************
	348	! Get variable id
	349	CALL check_err( nf90_inq_varid(ncid, "b", varid) )
	350	! Get the variable
	351	CALL check_err( nf90_get_var(ncid, varid, pt_b) )
[1150]	352
[1179]	353	! ++) Read p0 : reference pressure
	354	!****************************************************************************************
	355	! Get variable id
	356	CALL check_err( nf90_inq_varid(ncid, "p0", varid) )
	357	! Get the variable
	358	CALL check_err( nf90_get_var(ncid, varid, p0) )
[1150]	359
[524]	360
[1179]	361	ELSE ! old file
[524]	362
[1179]	363	! ++) Read the aerosol concentration month by month and concatenate to total variable varyear
	364	!****************************************************************************************
[1270]	365	DO imth=1, 12
[1179]	366	IF (imth.EQ.1) THEN
	367	cvar=TRIM(varname)//'JAN'
	368	ELSE IF (imth.EQ.2) THEN
	369	cvar=TRIM(varname)//'FEB'
	370	ELSE IF (imth.EQ.3) THEN
	371	cvar=TRIM(varname)//'MAR'
	372	ELSE IF (imth.EQ.4) THEN
	373	cvar=TRIM(varname)//'APR'
	374	ELSE IF (imth.EQ.5) THEN
	375	cvar=TRIM(varname)//'MAY'
	376	ELSE IF (imth.EQ.6) THEN
	377	cvar=TRIM(varname)//'JUN'
	378	ELSE IF (imth.EQ.7) THEN
	379	cvar=TRIM(varname)//'JUL'
	380	ELSE IF (imth.EQ.8) THEN
	381	cvar=TRIM(varname)//'AUG'
	382	ELSE IF (imth.EQ.9) THEN
	383	cvar=TRIM(varname)//'SEP'
	384	ELSE IF (imth.EQ.10) THEN
	385	cvar=TRIM(varname)//'OCT'
	386	ELSE IF (imth.EQ.11) THEN
	387	cvar=TRIM(varname)//'NOV'
	388	ELSE IF (imth.EQ.12) THEN
	389	cvar=TRIM(varname)//'DEC'
	390	END IF
	391
	392	! Get variable id
	393	CALL check_err( nf90_inq_varid(ncid, TRIM(cvar), varid) )
	394
	395	! Get the variable
	396	CALL check_err( nf90_get_var(ncid, varid, varmth) )
	397
	398	! Store in variable for the whole year
	399	varyear(:,:,:,imth)=varmth(:,:,:)
	400
	401	END DO
[1150]	402
[1179]	403	! Putting dummy
	404	psurf_glo2D(:,:,:) = not_valid
	405	load_glo2D(:,:,:) = not_valid
	406	pt_ap(:) = not_valid
	407	pt_b(:) = not_valid
[524]	408
[1179]	409	END IF
[524]	410
[1179]	411	! 4) Close file
	412	!****************************************************************************************
	413	CALL check_err( nf90_close(ncid) )
	414
[524]	415
[1179]	416	! 5) Transform the global field from 2D(iim, jjp+1) to 1D(klon_glo)
	417	!****************************************************************************************
	418	! Test if vertical levels have to be inversed
[782]	419
[1179]	420	IF ((pt_b(1) < pt_b(klev_src)) .OR. .NOT. new_file) THEN
[1183]	421	WRITE(lunout,*) 'Vertical axis in file ',TRIM(fname), ' needs to be inverted'
[1179]	422	WRITE(lunout,*) 'before pt_ap = ', pt_ap
	423	WRITE(lunout,*) 'before pt_b = ', pt_b
	424
	425	! Inverse vertical levels for varyear
[1270]	426	DO imth=1, 12
[1179]	427	varmth(:,:,:) = varyear(:,:,:,imth) ! use varmth temporarly
	428	DO k=1, klev_src
	429	DO j=1, jjm+1
	430	DO i=1,iim
	431	varyear(i,j,k,imth) = varmth(i,j,klev_src+1-k)
	432	END DO
	433	END DO
	434	END DO
	435	END DO
	436
	437	! Inverte vertical axes for pt_ap and pt_b
	438	varktmp(:) = pt_ap(:)
	439	DO k=1, klev_src
	440	pt_ap(k) = varktmp(klev_src+1-k)
	441	END DO
[782]	442
[1179]	443	varktmp(:) = pt_b(:)
	444	DO k=1, klev_src
	445	pt_b(k) = varktmp(klev_src+1-k)
	446	END DO
	447	WRITE(lunout,*) 'after pt_ap = ', pt_ap
	448	WRITE(lunout,*) 'after pt_b = ', pt_b
[524]	449
[1179]	450	ELSE
	451	WRITE(lunout,*) 'Vertical axis in file ',TRIM(fname), ' is ok, no vertical inversion is done'
	452	WRITE(lunout,*) 'pt_ap = ', pt_ap
	453	WRITE(lunout,*) 'pt_b = ', pt_b
	454	END IF
[524]	455
[1183]	456	! - Invert latitudes if necessary
[1270]	457	DO imth=1, 12
[1183]	458	IF (invert_lat) THEN
	459
	460	! Invert latitudes for the variable
	461	varmth(:,:,:) = varyear(:,:,:,imth) ! use varmth temporarly
	462	DO k=1,klev_src
	463	DO j=1,jjm+1
	464	DO i=1,iim
	465	varyear(i,j,k,imth) = varmth(i,jjm+1+1-j,k)
	466	END DO
	467	END DO
	468	END DO
	469
	470	! Invert latitudes for surface pressure
	471	vartmp(:,:) = psurf_glo2D(:,:,imth)
	472	DO j=1, jjm+1
[1179]	473	DO i=1,iim
[1183]	474	psurf_glo2D(i,j,imth)= vartmp(i,jjm+1+1-j)
[1179]	475	END DO
	476	END DO
[1183]	477
	478	! Invert latitudes for the load
	479	vartmp(:,:) = load_glo2D(:,:,imth)
	480	DO j=1, jjm+1
	481	DO i=1,iim
	482	load_glo2D(i,j,imth)= vartmp(i,jjm+1+1-j)
	483	END DO
[1179]	484	END DO
[1183]	485	END IF ! invert_lat
	486
[1179]	487	! Do zonal mead at poles and distribut at whole first and last latitude
	488	DO k=1, klev_src
	489	npole=0. ! North pole, j=1
	490	spole=0. ! South pole, j=jjm+1
	491	DO i=1,iim
	492	npole = npole + varyear(i,1,k,imth)
	493	spole = spole + varyear(i,jjm+1,k,imth)
	494	END DO
[1299]	495	npole = npole/REAL(iim)
	496	spole = spole/REAL(iim)
[1179]	497	varyear(:,1, k,imth) = npole
	498	varyear(:,jjm+1,k,imth) = spole
	499	END DO
	500	END DO ! imth
[1183]	501
[1270]	502	ALLOCATE(varyear_glo1D(klon_glo, klev_src, 12), stat=ierr)
[1179]	503	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 3',1)
	504
	505	! Transform from 2D to 1D field
	506	CALL grid2Dto1D_glo(varyear,varyear_glo1D)
	507	CALL grid2Dto1D_glo(psurf_glo2D,psurf_glo1D)
	508	CALL grid2Dto1D_glo(load_glo2D,load_glo1D)
[1249]	509
	510	ELSE
	511	ALLOCATE(varyear_glo1D(0,0,0))
	512	END IF ! is_mpi_root .AND. is_omp_root
	513
[1189]	514	!$OMP BARRIER
[1179]	515
	516	! 6) Distribute to all processes
	517	! Scatter global field(klon_glo) to local process domain(klon)
	518	! and distribute klev_src to all processes
	519	!****************************************************************************************
[524]	520
[1179]	521	! Distribute klev_src
	522	CALL bcast(klev_src)
[524]	523
[1179]	524	! Allocate and distribute pt_ap and pt_b
	525	IF (.NOT. ASSOCIATED(pt_ap)) THEN ! if pt_ap is allocated also pt_b is allocated
	526	ALLOCATE(pt_ap(klev_src), pt_b(klev_src), stat=ierr)
	527	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 4',1)
	528	END IF
	529	CALL bcast(pt_ap)
	530	CALL bcast(pt_b)
[524]	531
[1179]	532	! Allocate space for output pointer variable at local process
	533	IF (ASSOCIATED(pt_year)) DEALLOCATE(pt_year)
[1270]	534	ALLOCATE(pt_year(klon, klev_src, 12), stat=ierr)
	535	IF (ierr /= 0) CALL abort_gcm('get_aero_fromfile', 'pb in allocation 5',1)
[524]	536
[1179]	537	! Scatter global field to local domain at local process
	538	CALL scatter(varyear_glo1D, pt_year)
[1270]	539	CALL scatter(psurf_glo1D, psurf_out)
	540	CALL scatter(load_glo1D, load_out)
[524]	541
[1179]	542	! 7) Test for negative values
	543	!****************************************************************************************
	544	IF (MINVAL(pt_year) < 0.) THEN
	545	WRITE(lunout,*) 'Warning! Negative values read from file :', fname
	546	END IF
[524]	547
[1179]	548	END SUBROUTINE get_aero_fromfile
[524]	549
	550
[1179]	551	SUBROUTINE check_err(status)
	552	USE netcdf
	553	IMPLICIT NONE
[524]	554
[1179]	555	INCLUDE "iniprint.h"
	556	INTEGER, INTENT (IN) :: status
[524]	557
[1179]	558	IF (status /= NF90_NOERR) THEN
	559	WRITE(lunout,*) 'Error in get_aero_fromfile ',status
	560	CALL abort_gcm('get_aero_fromfile',trim(nf90_strerror(status)),1)
	561	END IF
[524]	562
[1179]	563	END SUBROUTINE check_err
	564
	565
	566	END MODULE readaerosol_mod

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

Download in other formats:

Original Format