Context Navigation

interp_sulf_input.F90 @ 5160

Last change on this file since 5160 was 3798, checked in by lguez, 4 years ago
Sync latest trunk changes to Ocean_skin
Property svn:keywords set to `Id`
File size: 9.5 KB

Line
1	!
2	! $Id: interp_sulf_input.F90 3798 2021-01-11 22:24:08Z abarral $
3	!
4	SUBROUTINE interp_sulf_input(debutphy,pdtphys,paprs,tr_seri)
5
6	USE netcdf95, ONLY: nf95_close, nf95_gw_var, nf95_inq_dimid, &
7	nf95_inq_varid, nf95_inquire_dimension, nf95_open
8	USE netcdf, ONLY: nf90_get_var, nf90_noerr, nf90_nowrite
9
10	USE mod_grid_phy_lmdz
11	USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root
12	USE mod_phys_lmdz_omp_data, ONLY : is_omp_root
13	USE phys_local_var_mod, ONLY : budg_3D_backgr_ocs, budg_3D_backgr_so2
14	USE phys_local_var_mod, ONLY : OCS_lifetime, SO2_lifetime
15	USE mod_phys_lmdz_para
16	USE dimphy
17	USE phys_cal_mod
18	USE infotrac_phy
19	USE aerophys
20	USE YOMCST
21
22	IMPLICIT NONE
23
24	include "dimensions.h"
25
26	! Variable input
27	REAL paprs(klon,klev+1)
28	REAL tr_seri(klon,klev,nbtr)
29	REAL, INTENT(IN) :: pdtphys ! Pas d'integration pour la physique (seconde)
30	LOGICAL, INTENT(IN) :: debutphy ! le flag de l'initialisation de la physique
31
32	! Variables locales
33	INTEGER n_lat ! number of latitudes in the input data
34	INTEGER n_lon ! number of longitudes in the input data
35	INTEGER, SAVE :: n_lev ! number of levels in the input data
36	!$OMP THREADPRIVATE(n_lev)
37	INTEGER n_mth ! number of months in the input data
38	INTEGER, SAVE :: mth_pre
39	!$OMP THREADPRIVATE(mth_pre)
40
41	! Champs reconstitues
42	REAL paprs_glo(klon_glo,klev+1)
43
44	REAL, POINTER:: latitude(:)
45	! (of input data sorted in strictly ascending order)
46
47	REAL, POINTER:: longitude(:)
48	! (of input data sorted in strictly ascending order)
49
50	REAL, POINTER:: time(:)
51	! (of input data sorted in strictly ascending order)
52
53	REAL, POINTER:: lev(:)
54	! levels of input data
55
56	REAL, ALLOCATABLE :: OCS_clim_in(:, :, :, :)
57	REAL, ALLOCATABLE :: SO2_clim_in(:, :, :, :)
58	REAL, ALLOCATABLE :: OCS_clim_mth(:, :, :)
59	REAL, ALLOCATABLE :: SO2_clim_mth(:, :, :)
60	REAL, ALLOCATABLE :: OCS_clim_tmp(:, :)
61	REAL, ALLOCATABLE :: SO2_clim_tmp(:, :)
62	REAL OCS_clim_glo(klon_glo,klev)
63	REAL SO2_clim_glo(klon_glo,klev)
64	REAL, ALLOCATABLE :: OCS_lifetime_in(:, :, :, :)
65	REAL, ALLOCATABLE :: SO2_lifetime_in(:, :, :, :)
66	REAL, ALLOCATABLE :: OCS_lifetime_mth(:, :, :)
67	REAL, ALLOCATABLE :: SO2_lifetime_mth(:, :, :)
68	REAL, ALLOCATABLE :: OCS_lifetime_tmp(:, :)
69	REAL, ALLOCATABLE :: SO2_lifetime_tmp(:, :)
70	REAL OCS_lifetime_glo(klon_glo,klev)
71	REAL SO2_lifetime_glo(klon_glo,klev)
72	!
73	REAL, ALLOCATABLE, SAVE :: OCS_clim(:,:)
74	REAL, ALLOCATABLE, SAVE :: SO2_clim(:,:)
75	!$OMP THREADPRIVATE(OCS_clim,SO2_clim)
76	!
77	INTEGER i, k, kk, j
78	REAL p_bound
79
80	! For NetCDF:
81	INTEGER ncid_in ! IDs for input files
82	INTEGER varid, ncerr
83
84	INTEGER, PARAMETER :: lev_input=17
85	!--pressure at interfaces of input data (in Pa)
86	REAL, DIMENSION(lev_input+1), PARAMETER :: &
87	paprs_input=(/ &
88	1.00000002e+05, 6.06530673e+04, 3.67879449e+04, &
89	2.23130165e+04, 1.35335286e+04, 8.20850004e+03, &
90	4.97870695e+03, 3.01973841e+03, 1.83156393e+03, &
91	1.11089968e+03, 6.73794715e+02, 4.08677153e+02, &
92	2.47875223e+02, 1.50343923e+02, 9.11881985e+01, &
93	5.53084382e+01, 3.35462635e+01, 0.0 /)
94	!
95	IF (.NOT.ALLOCATED(OCS_clim)) ALLOCATE(OCS_clim(klon,klev))
96	IF (.NOT.ALLOCATED(SO2_clim)) ALLOCATE(SO2_clim(klon,klev))
97
98	IF (debutphy.OR.mth_cur.NE.mth_pre) THEN
99
100	!--preparation of global fields
101	CALL gather(paprs, paprs_glo)
102
103	IF (is_mpi_root.AND.is_omp_root) THEN
104
105	!--reading emission files
106	CALL nf95_open("ocs_so2_annual_lmdz.nc", nf90_nowrite, ncid_in)
107
108	CALL nf95_inq_varid(ncid_in, "LEV", varid)
109	CALL nf95_gw_var(ncid_in, varid, lev)
110	n_lev = size(lev)
111
112	CALL nf95_inq_varid(ncid_in, "lat", varid)
113	CALL nf95_gw_var(ncid_in, varid, latitude)
114	n_lat = size(latitude)
115
116	CALL nf95_inq_varid(ncid_in, "lon", varid)
117	CALL nf95_gw_var(ncid_in, varid, longitude)
118	n_lon = size(longitude)
119
120	CALL nf95_inq_varid(ncid_in, "TIME", varid)
121	CALL nf95_gw_var(ncid_in, varid, time)
122	n_mth = size(time)
123
124	IF (.NOT.ALLOCATED(OCS_clim_in)) ALLOCATE(OCS_clim_in(n_lon, n_lat, n_lev, n_mth))
125	IF (.NOT.ALLOCATED(SO2_clim_in)) ALLOCATE(SO2_clim_in(n_lon, n_lat, n_lev, n_mth))
126	IF (.NOT.ALLOCATED(OCS_lifetime_in)) ALLOCATE(OCS_lifetime_in(n_lon, n_lat, n_lev, n_mth))
127	IF (.NOT.ALLOCATED(SO2_lifetime_in)) ALLOCATE(SO2_lifetime_in(n_lon, n_lat, n_lev, n_mth))
128
129	CALL nf95_inq_varid(ncid_in, "OCS", varid)
130	ncerr = nf90_get_var(ncid_in, varid, OCS_clim_in)
131	print *,'code erreur OCS=', ncerr, varid
132
133	CALL nf95_inq_varid(ncid_in, "SO2", varid)
134	ncerr = nf90_get_var(ncid_in, varid, SO2_clim_in)
135	print *,'code erreur SO2=', ncerr, varid
136
137	CALL nf95_inq_varid(ncid_in, "OCS_LIFET", varid)
138	ncerr = nf90_get_var(ncid_in, varid, OCS_lifetime_in)
139	print *,'code erreur OCS_lifetime_in=', ncerr, varid
140
141	CALL nf95_inq_varid(ncid_in, "SO2_LIFET", varid)
142	ncerr = nf90_get_var(ncid_in, varid, SO2_lifetime_in)
143	print *,'code erreur SO2_lifetime_in=', ncerr, varid
144
145	CALL nf95_close(ncid_in)
146
147	IF (.NOT.ALLOCATED(OCS_clim_mth)) ALLOCATE(OCS_clim_mth(n_lon, n_lat, n_lev))
148	IF (.NOT.ALLOCATED(SO2_clim_mth)) ALLOCATE(SO2_clim_mth(n_lon, n_lat, n_lev))
149	IF (.NOT.ALLOCATED(OCS_clim_tmp)) ALLOCATE(OCS_clim_tmp(klon_glo, n_lev))
150	IF (.NOT.ALLOCATED(SO2_clim_tmp)) ALLOCATE(SO2_clim_tmp(klon_glo, n_lev))
151	IF (.NOT.ALLOCATED(OCS_lifetime_mth)) ALLOCATE(OCS_lifetime_mth(n_lon, n_lat, n_lev))
152	IF (.NOT.ALLOCATED(SO2_lifetime_mth)) ALLOCATE(SO2_lifetime_mth(n_lon, n_lat, n_lev))
153	IF (.NOT.ALLOCATED(OCS_lifetime_tmp)) ALLOCATE(OCS_lifetime_tmp(klon_glo, n_lev))
154	IF (.NOT.ALLOCATED(SO2_lifetime_tmp)) ALLOCATE(SO2_lifetime_tmp(klon_glo, n_lev))
155
156	!---select the correct month, undo multiplication with 1.e12 (precision reasons)
157	!---correct latitudinal order and convert input from volume mixing ratio to mass mixing ratio
158	DO j=1,n_lat
159	SO2_clim_mth(:,j,:) = 1.e-12SO2_clim_in(:,n_lat+1-j,:,mth_cur)mSO2mol/mAIRmol
160	OCS_clim_mth(:,j,:) = 1.e-12OCS_clim_in(:,n_lat+1-j,:,mth_cur)mOCSmol/mAIRmol
161	SO2_lifetime_mth(:,j,:) = SO2_lifetime_in(:,n_lat+1-j,:,mth_cur)
162	OCS_lifetime_mth(:,j,:) = OCS_lifetime_in(:,n_lat+1-j,:,mth_cur)
163	ENDDO
164
165	!---reduce to a klon_glo grid but keep the levels
166	CALL grid2dTo1d_glo(OCS_clim_mth,OCS_clim_tmp)
167	CALL grid2dTo1d_glo(SO2_clim_mth,SO2_clim_tmp)
168	CALL grid2dTo1d_glo(OCS_lifetime_mth,OCS_lifetime_tmp)
169	CALL grid2dTo1d_glo(SO2_lifetime_mth,SO2_lifetime_tmp)
170
171	!--set lifetime to very high value in uninsolated areas
172	DO i=1, klon_glo
173	DO kk=1, n_lev
174	IF (OCS_lifetime_tmp(i,kk)==0.0) THEN
175	OCS_lifetime_tmp(i,kk)=1.0e12
176	ENDIF
177	IF (SO2_lifetime_tmp(i,kk)==0.0) THEN
178	SO2_lifetime_tmp(i,kk)=1.0e12
179	ENDIF
180	ENDDO
181	ENDDO
182
183	!---regrid weighted lifetime and climatologies
184	DO i=1, klon_glo
185	DO k=1, klev
186	OCS_lifetime_glo(i,k)=0.0
187	SO2_lifetime_glo(i,k)=0.0
188	OCS_clim_glo(i,k)=0.0
189	SO2_clim_glo(i,k)=0.0
190	DO kk=1, n_lev
191	OCS_lifetime_glo(i,k)=OCS_lifetime_glo(i,k)+ &
192	MAX(0.0,MIN(paprs_glo(i,k),paprs_input(kk))-MAX(paprs_glo(i,k+1),paprs_input(kk+1))) &
193	*OCS_lifetime_tmp(i,kk)/(paprs_glo(i,k)-paprs_glo(i,k+1))
194	SO2_lifetime_glo(i,k)=SO2_lifetime_glo(i,k)+ &
195	MAX(0.0,MIN(paprs_glo(i,k),paprs_input(kk))-MAX(paprs_glo(i,k+1),paprs_input(kk+1))) &
196	*SO2_lifetime_tmp(i,kk)/(paprs_glo(i,k)-paprs_glo(i,k+1))
197	OCS_clim_glo(i,k)=OCS_clim_glo(i,k)+ &
198	MAX(0.0,MIN(paprs_glo(i,k),paprs_input(kk))-MAX(paprs_glo(i,k+1),paprs_input(kk+1))) &
199	*OCS_clim_tmp(i,kk)/(paprs_glo(i,k)-paprs_glo(i,k+1))
200	SO2_clim_glo(i,k)=SO2_clim_glo(i,k)+ &
201	MAX(0.0,MIN(paprs_glo(i,k),paprs_input(kk))-MAX(paprs_glo(i,k+1),paprs_input(kk+1))) &
202	*SO2_clim_tmp(i,kk)/(paprs_glo(i,k)-paprs_glo(i,k+1))
203	ENDDO
204	ENDDO
205	ENDDO
206
207	ENDIF !--is_mpi_root and is_omp_root
208
209	!--keep memory of previous month
210	mth_pre=mth_cur
211
212	!--scatter global fields around
213	CALL scatter(OCS_clim_glo, OCS_clim)
214	CALL scatter(SO2_clim_glo, SO2_clim)
215	CALL scatter(OCS_lifetime_glo, OCS_lifetime)
216	CALL scatter(SO2_lifetime_glo, SO2_lifetime)
217
218	IF (is_mpi_root.AND.is_omp_root) THEN
219	!
220	DEALLOCATE(OCS_clim_in,SO2_clim_in)
221	DEALLOCATE(OCS_clim_mth,SO2_clim_mth)
222	DEALLOCATE(OCS_clim_tmp,SO2_clim_tmp)
223	DEALLOCATE(OCS_lifetime_in,SO2_lifetime_in)
224	DEALLOCATE(OCS_lifetime_mth,SO2_lifetime_mth)
225	DEALLOCATE(OCS_lifetime_tmp,SO2_lifetime_tmp)
226	!
227	ENDIF !--is_mpi_root and is_omp_root
228
229	ENDIF ! debutphy.OR.new month
230
231	!--set to background value everywhere in the very beginning, later only in the troposphere
232	!--a little dangerous as the MAXVAL is not computed on the global field
233	IF (debutphy.AND.MAXVAL(tr_seri).LT.1.e-30) THEN
234	p_bound=0.0
235	ELSE
236	p_bound=50000.
237	ENDIF
238
239	!--regridding tracer concentration on the vertical
240	DO i=1, klon
241	DO k=1, klev
242	!
243	!--OCS and SO2 prescribed back to their clim values below p_bound
244	IF (paprs(i,k).GT.p_bound) THEN
245	budg_3D_backgr_ocs(i,k)=OCS_clim(i,k)-tr_seri(i,k,id_OCS_strat)
246	budg_3D_backgr_so2(i,k)=SO2_clim(i,k)-tr_seri(i,k,id_SO2_strat)
247	tr_seri(i,k,id_OCS_strat)=OCS_clim(i,k)
248	tr_seri(i,k,id_SO2_strat)=SO2_clim(i,k)
249	ENDIF
250	ENDDO
251	ENDDO
252
253	!convert SO2_backgr_tend from kg(SO2)/kgA to kg(S)/m2/layer/s for saving as diagnostic
254	DO i=1, klon
255	DO k=1, klev
256	budg_3D_backgr_ocs(i,k)=budg_3D_backgr_ocs(i,k)mSatom/mOCSmol(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
257	budg_3D_backgr_so2(i,k)=budg_3D_backgr_so2(i,k)mSatom/mSO2mol(paprs(i,k)-paprs(i,k+1))/RG/pdtphys
258	ENDDO
259	ENDDO
260
261	RETURN
262
263	END SUBROUTINE interp_sulf_input

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

Download in other formats:

Original Format