Context Navigation

limit_netcdf.F90 @ 1323

Last change on this file since 1323 was 1323, checked in by Laurent Fairhead, 14 years ago

Changes made in r1293 are integrated into the trunk
Start files are identical between r1293 and this version

Les modifications de la r1293 sont intégrées à la trunk
Les fichiers start et startphy sont identiques entre la version 1293 et celle-ci

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

File size: 25.6 KB

Line
1	!
2	! $Id: limit_netcdf.F90 1323 2010-03-12 16:19:12Z fairhead $
3	!-------------------------------------------------------------------------------
4	!
5	SUBROUTINE limit_netcdf(interbar, extrap, oldice, masque)
6	!
7	!-------------------------------------------------------------------------------
8	! Author : L. Fairhead, 27/01/94
9	!-------------------------------------------------------------------------------
10	! Purpose: Boundary conditions files building for new model using climatologies.
11	! Both grids have to be regular.
12	!-------------------------------------------------------------------------------
13	! Note: This routine is designed to work for Earth
14	!-------------------------------------------------------------------------------
15	! Modification history:
16	! * 23/03/1994: Z. X. Li
17	! * 09/1999: L. Fairhead (netcdf reading in LMDZ.3.3)
18	! * 07/2001: P. Le Van
19	! * 11/2009: L. Guez (ozone day & night climatos, see etat0_netcdf.F90)
20	! * 12/2009: D. Cugnet (f77->f90, calendars, files from coupled runs)
21	!-------------------------------------------------------------------------------
22	#ifdef CPP_EARTH
23	USE dimphy
24	USE ioipsl, ONLY : ioget_year_len
25	USE phys_state_var_mod, ONLY : pctsrf
26	USE netcdf, ONLY : NF90_OPEN, NF90_CREATE, NF90_CLOSE, &
27	NF90_DEF_DIM, NF90_DEF_VAR, NF90_PUT_VAR, NF90_PUT_ATT, &
28	NF90_NOERR, NF90_NOWRITE, NF90_DOUBLE, NF90_GLOBAL, &
29	NF90_CLOBBER, NF90_ENDDEF, NF90_UNLIMITED
30	USE inter_barxy_m, only: inter_barxy
31	#endif
32	IMPLICIT NONE
33	!-------------------------------------------------------------------------------
34	! Arguments:
35	#include "dimensions.h"
36	#include "paramet.h"
37	#include "iniprint.h"
38	LOGICAL, INTENT(IN) :: interbar ! barycentric interpolation
39	LOGICAL, INTENT(IN) :: extrap ! SST extrapolation flag
40	LOGICAL, INTENT(IN) :: oldice ! old way ice computation
41	REAL, DIMENSION(iip1,jjp1), INTENT(IN) :: masque ! land mask
42	#ifndef CPP_EARTH
43	WRITE(lunout,*)'limit_netcdf: Earth-specific routine, needs Earth physics'
44	#else
45	!-------------------------------------------------------------------------------
46	! Local variables:
47	#include "control.h"
48	#include "logic.h"
49	#include "comvert.h"
50	#include "comgeom2.h"
51	#include "comconst.h"
52	#include "indicesol.h"
53
54	!--- For fractionary sub-cell use (old coding used soil index: 0,1,2,3) --------
55	LOGICAL, PARAMETER :: fracterre=.TRUE.
56
57	!--- INPUT NETCDF FILES NAMES --------------------------------------------------
58	CHARACTER(LEN=25) :: icefile, sstfile, dumstr
59	CHARACTER(LEN=25), PARAMETER :: famipsst='amipbc_sst_1x1.nc ', &
60	famipsic='amipbc_sic_1x1.nc ', &
61	fclimsst='amipbc_sst_1x1_clim.nc ', &
62	fclimsic='amipbc_sic_1x1_clim.nc ', &
63	fcpldsst='cpl_atm_sst.nc ', &
64	fcpldsic='cpl_atm_sic.nc ', &
65	frugo ='Rugos.nc ', &
66	falbe ='Albedo.nc '
67
68	!--- OUTPUT VARIABLES FOR NETCDF FILE ------------------------------------------
69	REAL, DIMENSION(klon) :: fi_ice, verif
70	REAL, DIMENSION(:,:), POINTER :: phy_rug=>NULL(), phy_ice=>NULL()
71	REAL, DIMENSION(:,:), POINTER :: phy_sst=>NULL(), phy_alb=>NULL()
72	REAL, DIMENSION(:,:), ALLOCATABLE :: phy_bil
73	REAL, DIMENSION(:,:,:), ALLOCATABLE :: pctsrf_t
74	INTEGER :: nbad
75
76	!--- VARIABLES FOR OUTPUT FILE WRITING -----------------------------------------
77	INTEGER :: ierr, nid, ndim, ntim, k
78	INTEGER, DIMENSION(2) :: dims
79	INTEGER :: id_tim, id_SST, id_BILS, id_RUG, id_ALB
80	INTEGER :: id_FOCE, id_FSIC, id_FTER, id_FLIC
81	INTEGER :: NF90_FORMAT
82	LOGICAL :: lCPL !--- T: IPCC-IPSL cpl model output files
83
84	!--- MICS (PHYSICAL KEYS, TIME) ------------------------------------------------
85	! INTEGER, PARAMETER :: longcles=20
86	! REAL, DIMENSION(longcles) :: clesphy0
87	INTEGER :: ndays
88
89	!--- INITIALIZATIONS -----------------------------------------------------------
90	#ifdef NC_DOUBLE
91	NF90_FORMAT=NF90_DOUBLE
92	#else
93	NF90_FORMAT=NF90_FLOAT
94	#endif
95
96	! CALL conf_gcm(99, .TRUE., clesphy0)
97	pi = 4.*ATAN(1.)
98	rad = 6371229.
99	daysec= 86400.
100	omeg = 2.*pi/daysec
101	g = 9.8
102	kappa = 0.2857143
103	cpp = 1004.70885
104	dtvr = daysec/FLOAT(day_step)
105	CALL inigeom
106
107	!--- Beware: anneeref (from gcm.def) is used to determine output time sampling
108	ndays=ioget_year_len(anneeref)
109
110	!--- RUGOSITY TREATMENT --------------------------------------------------------
111	WRITE(lunout,*) 'Traitement de la rugosite'
112	CALL get_2Dfield(frugo,'RUG',interbar,ndays,phy_rug,mask=masque(1:iim,:))
113
114	!--- OCEAN TREATMENT -----------------------------------------------------------
115	PRINT*, 'Traitement de la glace oceanique' ; icefile=''; lCPL=.FALSE.
116
117	! Input SIC file selection
118	icefile='fake'
119	IF(NF90_OPEN(famipsic,NF90_NOWRITE,nid)==NF90_NOERR) icefile=TRIM(famipsic)
120	IF(NF90_OPEN(fclimsic,NF90_NOWRITE,nid)==NF90_NOERR) icefile=TRIM(fclimsic)
121	IF(NF90_OPEN(fcpldsic,NF90_NOWRITE,nid)==NF90_NOERR) icefile=TRIM(fcpldsic)
122	SELECT CASE(icefile)
123	CASE(famipsic); dumstr='Amip.'
124	CASE(fclimsic); dumstr='Amip climatologique.'
125	CASE(fcpldsic); dumstr='de sortie du modele couplé IPSL/IPCC.';lCPL=.TRUE.
126	CASE('fake'); CALL abort_gcm('limit_netcdf','Fichier SIC non reconnu.',1)
127	END SELECT
128	ierr=NF90_CLOSE(nid)
129	WRITE(lunout,*)'Pour la glace de mer a ete choisi un fichier '//TRIM(dumstr)
130
131	CALL get_2Dfield(icefile,'SIC',interbar,ndays,phy_ice,flag=oldice)
132
133	ALLOCATE(pctsrf_t(klon,nbsrf,ndays))
134	DO k=1,ndays
135	fi_ice=phy_ice(:,k)
136	WHERE(fi_ice>=1.0 ) fi_ice=1.0
137	WHERE(fi_ice<EPSFRA) fi_ice=0.0
138	IF(fracterre) THEN
139	pctsrf_t(:,is_ter,k)=pctsrf(:,is_ter) ! land soil
140	pctsrf_t(:,is_lic,k)=pctsrf(:,is_lic) ! land ice
141	IF(lCPL) THEN ! SIC=pICE*(1-LIC-TER)
142	pctsrf_t(:,is_sic,k)=fi_ice*(1-pctsrf(:,is_lic)-pctsrf(:,is_ter))
143	ELSE ! SIC=pICE-LIC
144	pctsrf_t(:,is_sic,k)=fi_ice-pctsrf_t(:,is_lic,k)
145	END IF
146	WHERE(pctsrf_t(:,is_sic,k)<=0) pctsrf_t(:,is_sic,k)=0.
147	WHERE(1.0-zmasq<EPSFRA)
148	pctsrf_t(:,is_sic,k)=0.0
149	pctsrf_t(:,is_oce,k)=0.0
150	ELSEWHERE
151	WHERE(pctsrf_t(:,is_sic,k)>=1.0-zmasq)
152	pctsrf_t(:,is_sic,k)=1.0-zmasq
153	pctsrf_t(:,is_oce,k)=0.0
154	ELSEWHERE
155	pctsrf_t(:,is_oce,k)=1.0-zmasq-pctsrf_t(:,is_sic,k)
156	WHERE(pctsrf_t(:,is_oce,k)<EPSFRA)
157	pctsrf_t(:,is_oce,k)=0.0
158	pctsrf_t(:,is_sic,k)=1.0-zmasq
159	END WHERE
160	END WHERE
161	END WHERE
162	nbad=COUNT(pctsrf_t(:,is_oce,k)<0.0)
163	IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb point = ',nbad
164	nbad=COUNT(abs(sum(pctsrf_t(:,:,k),dim=2)-1.0)>EPSFRA)
165	IF(nbad>0) WRITE(lunout,*) 'pb sous surface pour nb points = ',nbad
166	ELSE
167	pctsrf_t(:,is_ter,k)=pctsrf(:,is_ter)
168	WHERE(NINT(pctsrf(:,is_ter))==1)
169	pctsrf_t(:,is_sic,k)=0.
170	pctsrf_t(:,is_oce,k)=0.
171	WHERE(fi_ice>=1.e-5)
172	pctsrf_t(:,is_lic,k)=fi_ice
173	pctsrf_t(:,is_ter,k)=pctsrf_t(:,is_ter,k)-pctsrf_t(:,is_lic,k)
174	ELSEWHERE
175	pctsrf_t(:,is_lic,k)=0.0
176	END WHERE
177	ELSEWHERE
178	pctsrf_t(:,is_lic,k) = 0.0
179	WHERE(fi_ice>=1.e-5)
180	pctsrf_t(:,is_ter,k)=0.0
181	pctsrf_t(:,is_sic,k)=fi_ice
182	pctsrf_t(:,is_oce,k)=1.0-pctsrf_t(:,is_sic,k)
183	ELSEWHERE
184	pctsrf_t(:,is_sic,k)=0.0
185	pctsrf_t(:,is_oce,k)=1.0
186	END WHERE
187	END WHERE
188	verif=sum(pctsrf_t(:,:,k),dim=2)
189	nbad=COUNT(verif<1.0-1.e-5.OR.verif>1.0+1.e-5)
190	IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb point = ',nbad
191	END IF
192	END DO
193	DEALLOCATE(phy_ice)
194
195	!--- SST TREATMENT -------------------------------------------------------------
196	WRITE(lunout,*) 'Traitement de la sst' ; sstfile=''; lCPL=.FALSE.
197
198	! Input SST file selection
199	sstfile='fake'
200	IF(NF90_OPEN(famipsst,NF90_NOWRITE,nid)==NF90_NOERR) sstfile=TRIM(famipsst)
201	IF(NF90_OPEN(fclimsst,NF90_NOWRITE,nid)==NF90_NOERR) sstfile=TRIM(fclimsst)
202	IF(NF90_OPEN(fcpldsst,NF90_NOWRITE,nid)==NF90_NOERR) sstfile=TRIM(fcpldsst)
203	SELECT CASE(icefile)
204	CASE(famipsic); dumstr='Amip.'
205	CASE(fclimsic); dumstr='Amip climatologique.'
206	CASE(fcpldsic); dumstr='de sortie du modele couplé IPSL/IPCC.';lCPL=.TRUE.
207	CASE('fake'); CALL abort_gcm('limit_netcdf','Fichier SST non reconnu',1)
208	END SELECT
209	ierr=NF90_CLOSE(nid)
210	WRITE(lunout,*)'Pour la temperature de mer a ete choisi un fichier '//TRIM(dumstr)
211
212	CALL get_2Dfield(trim(sstfile),'SST',interbar,ndays,phy_sst,flag=extrap)
213
214	!--- ALBEDO TREATMENT ----------------------------------------------------------
215	WRITE(lunout,*) 'Traitement de l albedo'
216	CALL get_2Dfield(falbe,'ALB',interbar,ndays,phy_alb)
217
218	!--- REFERENCE GROUND HEAT FLUX TREATMENT --------------------------------------
219	ALLOCATE(phy_bil(klon,ndays)); phy_bil=0.0
220
221	!--- OUTPUT FILE WRITING -------------------------------------------------------
222	WRITE(lunout,*) 'Ecriture du fichier limit'
223
224	!--- File creation
225	ierr=NF90_CREATE("limit.nc",NF90_CLOBBER,nid)
226	ierr=NF90_PUT_ATT(nid,NF90_GLOBAL,"title","Fichier conditions aux limites")
227
228	!--- Dimensions creation
229	ierr=NF90_DEF_DIM(nid,"points_physiques",klon,ndim)
230	ierr=NF90_DEF_DIM(nid,"time",NF90_UNLIMITED,ntim)
231
232	dims=(/ndim,ntim/)
233
234	!--- Variables creation
235	ierr=NF90_DEF_VAR(nid,"TEMPS",NF90_FORMAT,dims,id_tim)
236	ierr=NF90_DEF_VAR(nid,"FOCE", NF90_FORMAT,dims,id_FOCE)
237	ierr=NF90_DEF_VAR(nid,"FSIC", NF90_FORMAT,dims,id_FSIC)
238	ierr=NF90_DEF_VAR(nid,"FTER", NF90_FORMAT,dims,id_FTER)
239	ierr=NF90_DEF_VAR(nid,"FLIC", NF90_FORMAT,dims,id_FLIC)
240	ierr=NF90_DEF_VAR(nid,"SST", NF90_FORMAT,dims,id_SST)
241	ierr=NF90_DEF_VAR(nid,"BILS", NF90_FORMAT,dims,id_BILS)
242	ierr=NF90_DEF_VAR(nid,"ALB", NF90_FORMAT,dims,id_ALB)
243	ierr=NF90_DEF_VAR(nid,"RUG", NF90_FORMAT,dims,id_RUG)
244
245	!--- Attributes creation
246	ierr=NF90_PUT_ATT(nid,id_tim, "title","Jour dans l annee")
247	ierr=NF90_PUT_ATT(nid,id_FOCE,"title","Fraction ocean")
248	ierr=NF90_PUT_ATT(nid,id_FSIC,"title","Fraction glace de mer")
249	ierr=NF90_PUT_ATT(nid,id_FTER,"title","Fraction terre")
250	ierr=NF90_PUT_ATT(nid,id_FLIC,"title","Fraction land ice")
251	ierr=NF90_PUT_ATT(nid,id_SST ,"title","Temperature superficielle de la mer")
252	ierr=NF90_PUT_ATT(nid,id_BILS,"title","Reference flux de chaleur au sol")
253	ierr=NF90_PUT_ATT(nid,id_ALB, "title","Albedo a la surface")
254	ierr=NF90_PUT_ATT(nid,id_RUG, "title","Rugosite")
255
256	ierr=NF90_ENDDEF(nid)
257
258	!--- Variables saving
259	ierr=NF90_PUT_VAR(nid,id_tim,(/(DBLE(k),k=1,ndays)/))
260	ierr=NF90_PUT_VAR(nid,id_FOCE,pctsrf_t(:,is_oce,:),(/1,1/),(/klon,ndays/))
261	ierr=NF90_PUT_VAR(nid,id_FSIC,pctsrf_t(:,is_sic,:),(/1,1/),(/klon,ndays/))
262	ierr=NF90_PUT_VAR(nid,id_FTER,pctsrf_t(:,is_ter,:),(/1,1/),(/klon,ndays/))
263	ierr=NF90_PUT_VAR(nid,id_FLIC,pctsrf_t(:,is_lic,:),(/1,1/),(/klon,ndays/))
264	ierr=NF90_PUT_VAR(nid,id_SST ,phy_sst(:,:),(/1,1/),(/klon,ndays/))
265	ierr=NF90_PUT_VAR(nid,id_BILS,phy_bil(:,:),(/1,1/),(/klon,ndays/))
266	ierr=NF90_PUT_VAR(nid,id_ALB ,phy_alb(:,:),(/1,1/),(/klon,ndays/))
267	ierr=NF90_PUT_VAR(nid,id_RUG ,phy_rug(:,:),(/1,1/),(/klon,ndays/))
268
269	ierr=NF90_CLOSE(nid)
270
271	DEALLOCATE(pctsrf_t,phy_sst,phy_bil,phy_alb,phy_rug)
272	#endif
273	! of #ifdef CPP_EARTH
274	STOP
275
276
277	!===============================================================================
278	!
279	CONTAINS
280	!
281	!===============================================================================
282
283
284	!-------------------------------------------------------------------------------
285	!
286	SUBROUTINE get_2Dfield(fnam, mode, ibar, ndays, champo, flag, mask)
287	!
288	!-------------------------------------------------------------------------------
289	! Comments:
290	! There are two assumptions concerning the NetCDF files, that are satisfied
291	! with files that are conforming NC convention:
292	! 1) The last dimension of the variables used is the time record.
293	! 2) Dimensional variables have the same names as corresponding dimensions.
294	!-------------------------------------------------------------------------------
295	USE netcdf, ONLY : NF90_OPEN, NF90_INQ_VARID, NF90_INQUIRE_VARIABLE, &
296	NF90_CLOSE, NF90_INQ_DIMID, NF90_INQUIRE_DIMENSION, &
297	NF90_GET_VAR, NF90_GET_ATT
298	USE dimphy, ONLY : klon
299	USE phys_state_var_mod, ONLY : pctsrf
300	IMPLICIT NONE
301	#include "dimensions.h"
302	#include "paramet.h"
303	#include "comgeom2.h"
304	#include "control.h"
305	#include "indicesol.h"
306	#include "iniprint.h"
307	!-------------------------------------------------------------------------------
308	! Arguments:
309	CHARACTER(LEN=*), INTENT(IN) :: fnam ! NetCDF file name
310	CHARACTER(LEN=3), INTENT(IN) :: mode ! RUG, SIC, SST or ALB
311	LOGICAL, INTENT(IN) :: ibar ! interp on pressure levels
312	INTEGER, INTENT(IN) :: ndays ! current year number of days
313	REAL, POINTER, DIMENSION(:,:) :: champo ! output field = f(t)
314	LOGICAL, OPTIONAL, INTENT(IN) :: flag
315	! flag=T means: extrapolation (SST case) or old ice (SIC case)
316	REAL, OPTIONAL, DIMENSION(iim,jjp1), INTENT(IN) :: mask
317	!-------------------------------------------------------------------------------
318	! Local variables:
319	!--- NetCDF
320	INTEGER :: ierr, ncid, varid ! NetCDF identifiers
321	CHARACTER(LEN=30) :: dnam ! dimension name
322	CHARACTER(LEN=80) :: varname ! NetCDF variable name
323	!--- dimensions
324	INTEGER, DIMENSION(4) :: dids ! NetCDF dimensions identifiers
325	REAL, ALLOCATABLE, DIMENSION(:) :: dlon_ini ! initial longitudes vector
326	REAL, ALLOCATABLE, DIMENSION(:) :: dlat_ini ! initial latitudes vector
327	REAL, POINTER, DIMENSION(:) :: dlon, dlat ! reordered lon/lat vectors
328	!--- fields
329	INTEGER :: imdep, jmdep, lmdep ! dimensions of 'champ'
330	REAL, ALLOCATABLE, DIMENSION(:,:) :: champ ! wanted field on initial grid
331	REAL, ALLOCATABLE, DIMENSION(:) :: yder, timeyear
332	REAL, DIMENSION(iim,jjp1) :: champint ! interpolated field
333	REAL, ALLOCATABLE, DIMENSION(:,:,:) :: champtime
334	REAL, ALLOCATABLE, DIMENSION(:,:,:) :: champan
335	REAL :: time, by
336	!--- input files
337	CHARACTER(LEN=20) :: cal_in ! calendar
338	INTEGER :: ndays_in ! number of days
339	!--- misc
340	INTEGER :: i, j, k, l ! loop counters
341	REAL, ALLOCATABLE, DIMENSION(:,:) :: work ! used for extrapolation
342	CHARACTER(LEN=25) :: title ! for messages
343	LOGICAL :: extrp ! flag for extrapolation
344	REAL :: chmin, chmax
345	!-------------------------------------------------------------------------------
346	!---Variables depending on keyword 'mode' --------------------------------------
347	NULLIFY(champo)
348	SELECT CASE(mode)
349	CASE('RUG'); varname='RUGOS'; title='Rugosite'
350	CASE('SIC'); varname='sicbcs'; title='Sea-ice'
351	CASE('SST'); varname='tosbcs'; title='SST'
352	CASE('ALB'); varname='ALBEDO'; title='Albedo'
353	END SELECT
354	extrp=(flag.AND.mode=='SST')
355
356	!--- GETTING SOME DIMENSIONAL VARIABLES FROM FILE ------------------------------
357	ierr=NF90_OPEN(fnam,NF90_NOWRITE,ncid); CALL ncerr(ierr,fnam)
358	ierr=NF90_INQ_VARID(ncid,varname,varid); CALL ncerr(ierr,fnam)
359	ierr=NF90_INQUIRE_VARIABLE(ncid,varid,dimids=dids); CALL ncerr(ierr,fnam)
360
361	!--- Longitude
362	ierr=NF90_INQUIRE_DIMENSION(ncid,dids(1),name=dnam,len=imdep)
363	CALL ncerr(ierr,fnam); ALLOCATE(dlon_ini(imdep),dlon(imdep))
364	ierr=NF90_INQ_VARID(ncid,dnam,varid); CALL ncerr(ierr,fnam)
365	ierr=NF90_GET_VAR(ncid,varid,dlon_ini); CALL ncerr(ierr,fnam)
366	WRITE(lunout,*) 'variable ', dnam, 'dimension ', imdep
367
368	!--- Latitude
369	ierr=NF90_INQUIRE_DIMENSION(ncid,dids(2),name=dnam,len=jmdep)
370	CALL ncerr(ierr,fnam); ALLOCATE(dlat_ini(jmdep),dlat(jmdep))
371	ierr=NF90_INQ_VARID(ncid,dnam,varid); CALL ncerr(ierr,fnam)
372	ierr=NF90_GET_VAR(ncid,varid,dlat_ini); CALL ncerr(ierr,fnam)
373	WRITE(lunout,*) 'variable ', dnam, 'dimension ', jmdep
374
375	!--- Time (variable is not needed - it is rebuilt - but calendar is)
376	ierr=NF90_INQUIRE_DIMENSION(ncid,dids(3),name=dnam,len=lmdep)
377	CALL ncerr(ierr,fnam); ALLOCATE(timeyear(lmdep))
378	ierr=NF90_INQ_VARID(ncid,dnam,varid); CALL ncerr(ierr,fnam)
379	cal_in=' '
380	ierr=NF90_GET_ATT(ncid,varid,'calendar',cal_in)
381	IF(ierr/=NF90_NOERR) THEN
382	SELECT CASE(mode)
383	CASE('RUG','ALB'); cal_in='360d'
384	CASE('SIC','SST'); cal_in='gregorian'
385	END SELECT
386	WRITE(lunout,*)'ATTENTION: variable ''time'' sans attribut ''calendrier'' d&
387	&ans '//TRIM(fnam)//'. On choisit la valeur par defaut.'
388	END IF
389	WRITE(lunout,*) 'variable ', dnam, 'dimension ', lmdep, 'calendrier ', cal_in
390
391	!--- CONSTRUCTING THE INPUT TIME VECTOR FOR INTERPOLATION ----------------------
392	!--- Determining input file number of days, depending on calendar
393	ndays_in=year_len(anneeref,cal_in)
394
395	!--- Time vector reconstruction (time vector from file is not trusted)
396	!--- If input records are not monthly, time sampling has to be constant !
397	timeyear=mid_months(anneeref,cal_in,lmdep)
398	IF(lmdep/=12) WRITE(lunout,'(a,i3,a)')'Note: les fichiers de '//TRIM(mode) &
399	//' ne comportent pas 12, mais ',lmdep,' renregistrements.'
400
401	!--- GETTING THE FIELD AND INTERPOLATING IT ------------------------------------
402	ALLOCATE(champ(imdep,jmdep),champtime(iim,jjp1,lmdep))
403	IF(extrp) ALLOCATE(work(imdep,jmdep))
404
405	WRITE(lunout,*)
406	WRITE(lunout,'(a,i3,a)')'LECTURE ET INTERPOLATION HORIZ. DE ',lmdep,' CHAMPS.'
407	ierr=NF90_INQ_VARID(ncid,varname,varid); CALL ncerr(ierr,fnam)
408	DO l=1,lmdep
409	ierr=NF90_GET_VAR(ncid,varid,champ,(/1,1,l/),(/imdep,jmdep,1/))
410	CALL ncerr(ierr,fnam)
411	CALL conf_dat2d(title,imdep,jmdep,dlon_ini,dlat_ini,dlon,dlat,champ,ibar)
412
413	IF(extrp) CALL extrapol(champ,imdep,jmdep,999999.,.TRUE.,.TRUE.,2,work)
414
415	IF(ibar.AND..NOT.(mode=='SIC'.AND.flag)) THEN
416	IF(l==1) THEN
417	WRITE(lunout,*) &
418	'-------------------------------------------------------------------------'
419	WRITE(lunout,*) &
420	'$$$ Utilisation de l''interpolation barycentrique pour '//TRIM(title)//' $$$'
421	WRITE(lunout,*) &
422	'-------------------------------------------------------------------------'
423	END IF
424	IF(mode=='RUG') champ=LOG(champ)
425	CALL inter_barxy(dlon, dlat(:jmdep-1), champ, rlonu(:iim), &
426	rlatv, champint)
427	IF(mode=='RUG') THEN
428	champint=EXP(champint)
429	WHERE(NINT(mask)/=1) champint=0.001
430	END IF
431	ELSE
432	SELECT CASE(mode)
433	CASE('RUG'); CALL rugosite(imdep, jmdep, dlon, dlat, champ, &
434	iim, jjp1, rlonv, rlatu, champint, mask)
435	CASE('SIC'); CALL sea_ice (imdep, jmdep, dlon, dlat, champ, &
436	iim, jjp1, rlonv, rlatu, champint)
437	CASE('SST','ALB'); CALL grille_m(imdep, jmdep, dlon, dlat, champ, &
438	iim, jjp1, rlonv, rlatu, champint)
439	END SELECT
440	END IF
441	champtime(:,:,l)=champint
442	END DO
443	ierr=NF90_CLOSE(ncid); CALL ncerr(ierr,fnam)
444
445	DEALLOCATE(dlon_ini,dlat_ini,dlon,dlat,champ)
446	IF(extrp) DEALLOCATE(work)
447
448	!--- TIME INTERPOLATION --------------------------------------------------------
449	WRITE(lunout,*)
450	WRITE(lunout,*)'INTERPOLATION TEMPORELLE.'
451	WRITE(lunout,"(2x,' Vecteur temps en entree: ',10f6.1)") timeyear
452	WRITE(lunout,"(2x,' Vecteur temps en sortie de 0 a ',i3)") ndays
453	ALLOCATE(yder(lmdep),champan(iip1,jjp1,ndays))
454	DO j=1,jjp1
455	DO i=1,iim
456	CALL spline(timeyear,champtime(i,j,:),lmdep,1.e30,1.e30,yder)
457	DO k=1,ndays
458	time=FLOAT((k-1)*ndays_in)/FLOAT(ndays)
459	CALL splint(timeyear,champtime(i,j,:),yder,lmdep,time,by)
460	champan(i,j,k) = by
461	END DO
462	END DO
463	END DO
464	champan(iip1,:,:)=champan(1,:,:)
465	DEALLOCATE(yder,champtime,timeyear)
466
467	!--- Checking the result
468	DO j=1,jjp1
469	CALL minmax(iip1,champan(1,j,10),chmin,chmax)
470	WRITE(lunout,*)' '//TRIM(title)//' au temps 10 ',chmin,chmax,j
471	END DO
472
473	!--- SPECIAL FILTER FOR SST: SST>271.38 ----------------------------------------
474	IF(mode=='SST') THEN
475	WRITE(lunout,*) 'Filtrage de la SST: SST >= 271.38'
476	WHERE(champan<271.38) champan=271.38
477	END IF
478
479	!--- SPECIAL FILTER FOR SIC: 0.0<SIC<1.0 ---------------------------------------
480	IF(mode=='SIC') THEN
481	WRITE(lunout,*) 'Filtrage de la SIC: 0.0 < Sea-ice < 1.0'
482	champan(:,:,:)=champan(:,:,:)/100.
483	champan(iip1,:,:)=champan(1,:,:)
484	WHERE(champan>1.0) champan=1.0
485	WHERE(champan<0.0) champan=0.0
486	END IF
487
488	write(,)'coin1'
489	!--- DYNAMICAL TO PHYSICAL GRID ------------------------------------------------
490	ALLOCATE(champo(klon,ndays))
491	DO k=1,ndays
492	CALL gr_dyn_fi(1,iip1,jjp1,klon,champan(1,1,k),champo(1,k))
493	END DO
494	write(,)'coin2'
495	DEALLOCATE(champan)
496	write(,)'coin3'
497	END SUBROUTINE get_2Dfield
498	!
499	!-------------------------------------------------------------------------------
500
501
502
503	!-------------------------------------------------------------------------------
504	!
505	FUNCTION year_len(y,cal_in)
506	!
507	!-------------------------------------------------------------------------------
508	USE ioipsl, ONLY : ioget_calendar,ioconf_calendar,lock_calendar,ioget_year_len
509	IMPLICIT NONE
510	!-------------------------------------------------------------------------------
511	! Arguments:
512	INTEGER :: year_len
513	INTEGER, INTENT(IN) :: y
514	CHARACTER(LEN=*), INTENT(IN) :: cal_in
515	!-------------------------------------------------------------------------------
516	! Local variables:
517	CHARACTER(LEN=20) :: cal_out ! calendar (for outputs)
518	!-------------------------------------------------------------------------------
519	!--- Getting the input calendar to reset at the end of the function
520	CALL ioget_calendar(cal_out)
521
522	!--- Unlocking calendar and setting it to wanted one
523	CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in))
524
525	!--- Getting the number of days in this year
526	year_len=ioget_year_len(y)
527
528	!--- Back to original calendar
529	CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out))
530
531	END FUNCTION year_len
532	!
533	!-------------------------------------------------------------------------------
534
535
536	!-------------------------------------------------------------------------------
537	!
538	FUNCTION mid_months(y,cal_in,nm)
539	!
540	!-------------------------------------------------------------------------------
541	USE ioipsl, ONLY : ioget_calendar,ioconf_calendar,lock_calendar,ioget_mon_len
542	IMPLICIT NONE
543	!-------------------------------------------------------------------------------
544	! Arguments:
545	INTEGER, INTENT(IN) :: y ! year
546	CHARACTER(LEN=*), INTENT(IN) :: cal_in ! calendar
547	INTEGER, INTENT(IN) :: nm ! months/year number
548	REAL, DIMENSION(nm) :: mid_months ! mid-month times
549	!-------------------------------------------------------------------------------
550	! Local variables:
551	CHARACTER(LEN=99) :: mess ! error message
552	CHARACTER(LEN=20) :: cal_out ! calendar (for outputs)
553	INTEGER, DIMENSION(nm) :: mnth ! months lengths (days)
554	INTEGER :: m ! months counter
555	INTEGER :: nd ! number of days
556	!-------------------------------------------------------------------------------
557	nd=year_len(y,cal_in)
558
559	IF(nm==12) THEN
560
561	!--- Getting the input calendar to reset at the end of the function
562	CALL ioget_calendar(cal_out)
563
564	!--- Unlocking calendar and setting it to wanted one
565	CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in))
566
567	!--- Getting the length of each month
568	DO m=1,nm; mnth(m)=ioget_mon_len(y,m); END DO
569
570	!--- Back to original calendar
571	CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out))
572
573	ELSE IF(MODULO(nd,nm)/=0) THEN
574	WRITE(mess,'(a,i3,a,i3,a)')'Unconsistent calendar: ',nd,' days/year, but ',&
575	nm,' months/year. Months number should divide days number.'
576	CALL abort_gcm('mid_months',TRIM(mess),1)
577
578	ELSE
579	mnth=(/(m,m=1,nm,nd/nm)/)
580	END IF
581
582	!--- Mid-months times
583	mid_months(1)=0.5*FLOAT(mnth(1))
584	DO k=2,nm
585	mid_months(k)=mid_months(k-1)+0.5*FLOAT(mnth(k-1)+mnth(k))
586	END DO
587
588	END FUNCTION mid_months
589	!
590	!-------------------------------------------------------------------------------
591
592
593	!-------------------------------------------------------------------------------
594	!
595	SUBROUTINE ncerr(ncres,fnam)
596	!
597	!-------------------------------------------------------------------------------
598	! Purpose: NetCDF errors handling.
599	!-------------------------------------------------------------------------------
600	USE netcdf, ONLY : NF90_NOERR, NF90_STRERROR
601	IMPLICIT NONE
602	!-------------------------------------------------------------------------------
603	! Arguments:
604	INTEGER, INTENT(IN) :: ncres
605	CHARACTER(LEN=*), INTENT(IN) :: fnam
606	!-------------------------------------------------------------------------------
607	#include "iniprint.h"
608	IF(ncres/=NF90_NOERR) THEN
609	WRITE(lunout,*)'Problem with file '//TRIM(fnam)//' in routine limit_netcdf.'
610	CALL abort_gcm('limit_netcdf',NF90_STRERROR(ncres),1)
611	END IF
612
613	END SUBROUTINE ncerr
614	!
615	!-------------------------------------------------------------------------------
616
617
618	END SUBROUTINE limit_netcdf

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

Download in other formats:

Original Format