source: LMDZ6/branches/LMDZ_ECRad/libf/dyn3d_common/conf_dat_m.F90 @ 5434

Last change on this file since 5434 was 2293, checked in by dcugnet, 10 years ago

Initial states creation routines have been reorganized and simplified.
As far as possible, dynamics and physics related routines have been
separated.
Some routines have been converted to fortran 90 and repeated codes sections
have been "factorized".
Array/vector arguments have become implicit in some routines to avoid usage
of "dimensions.h" ; possible for routines with explicit interfaces and if
iim and jjm can be deduced from arguments sizes.

  • dynlonlat_phylonlat/ce0l.F90 calls now phylmd/etat0phys_netcdf.F90 and dyn3d/etat0dyn_netcdf.F90 that replace phylmd/etat0_netcdf.F90. start.nc and startphy.nc creations are now independant.
  • startvar.F90 has been suppressed ; corresponding operations have been simplified and embedded in etat0*_netcdf.F90 routines as internal procedures.
  • Routines converted to fortran 90 and "factorized":
    • dyn3d_common/conf_dat_m.F90 (replaces dyn3d_common/conf_dat2d.F

and dyn3d_common/conf_dat3d.F)

  • dyn3d/dynredem.F90 (replaces dyn3d/dynredem.F)
  • dyn3d/dynetat0.F90 (replaces dyn3d/dynetat0.F)
  • phylmd/grid_noro_m.F90 (replaces dyn3d_common/grid_noro.F)
  • dynlonlat_phylonlat/grid_atob_m.F90 (replaces dyn3d_common/grid_atob.F)
  • dyn3d_common/caldyn0.F90 (replaces dyn3d_common/caldyn0.F)
  • dyn3d_common/covcont.F90 (replaces dyn3d_common/covcont.F)
  • dyn3d_common/pression.F90 (replaces dyn3d_common/pression.F)
  • phylmd/phyredem.F90 and phylmd/limit_netcdf.F90 have been slightly factorized.

TO DO:

  • little fix needed in grid_noro_m.F90 ; untouched yet to ensure results are exactly the same as before. Unsmoothed orography is used to compute "zphi", but smoothed (should be unsmoothed) one is used at poles.
  • add the dyn3dmem versions of dynredem.F90 and dynetat0.F90 (dynredem_loc.F90 and dynetat0_loc.F90, untested yet).
  • test compilation in parallel mode for a single processor.
File size: 8.9 KB
Line 
1MODULE conf_dat_m
2!
3!*******************************************************************************
4
5  PRIVATE
6  PUBLIC :: conf_dat2d, conf_dat3d
7
8
9CONTAINS
10
11
12!-------------------------------------------------------------------------------
13!
14SUBROUTINE conf_dat2d(title, xd, yd, xf, yf, champd, interbar)
15!
16!-------------------------------------------------------------------------------
17! Author: P. Le Van
18!-------------------------------------------------------------------------------
19! Purpose: Configure the 2D data field "champd" so that:
20!     - Longitudes are in [ -pi     pi   ]
21!     - Latitudes  are in [  pi/2. -pi/2 ]
22!  * xd / yd are initial lon / lats.
23!  * xf / yf are output  lon / lats, possibly modified to satisfy configuration.
24!  * interbar is TRUE for barycentric interpolation.
25!-------------------------------------------------------------------------------
26  USE assert_eq_m, ONLY: assert_eq
27  IMPLICIT NONE
28!-------------------------------------------------------------------------------
29! Arguments:
30  CHARACTER(LEN=*), INTENT(IN)    :: title
31  REAL,             INTENT(IN)    :: xd(:), yd(:)  ! dim (lons) (lats)
32  REAL,             INTENT(INOUT) :: xf(:), yf(:)  ! dim (lons) (lats)
33  REAL,             INTENT(INOUT) :: champd(:,:)   ! dim (lons,lats)
34  LOGICAL,          INTENT(IN)    :: interbar
35!-------------------------------------------------------------------------------
36! Local variables:
37  CHARACTER(LEN=256) :: modname="conf_dat2d"
38  INTEGER :: i, j, ip180, ind, lons, lats
39  LOGICAL :: radlon, invlon ,radlat, invlat
40  REAL    :: pi, pis2, depi, rlatmin, rlatmax, oldxd1
41  REAL, ALLOCATABLE :: xtemp(:) , ytemp(:), champf(:,:)
42!-------------------------------------------------------------------------------
43  lons=assert_eq(SIZE(xd),SIZE(xf),SIZE(champd,1),TRIM(modname)//" lons")
44  lats=assert_eq(SIZE(yd),SIZE(yf),SIZE(champd,2),TRIM(modname)//" lats")
45  ALLOCATE(xtemp(lons),ytemp(lats)); xtemp(:)=xd(:); ytemp(:)=yd(:)
46  ALLOCATE(champf(lons,lats))
47  pi = 2. * ASIN(1.)
48  pis2 = pi/2.
49  depi = 2. * pi
50  radlon=.FALSE.; invlon=.FALSE.
51  IF     (xtemp(1)>=-pi-0.5.AND.xtemp(lons)<=  pi+0.5) THEN
52    radlon = .TRUE.;  invlon = .FALSE.
53  ELSE IF(xtemp(1)>=   -0.5.AND.xtemp(lons)<=depi+0.5) THEN
54    radlon = .TRUE.;  invlon = .TRUE.
55  ELSE IF(xtemp(1)>= -180.5.AND.xtemp(lons)<=   180.5) THEN
56    radlon = .FALSE.; invlon = .FALSE.
57  ELSE IF(xtemp(1)>=   -0.5.AND.xtemp(lons)<=   360.5) THEN
58    radlon = .FALSE.; invlon = .TRUE.
59  ELSE; WRITE(6,*) 'Problems with data longitudes for file '//TRIM(title)
60  END IF
61  invlat = ytemp(1)<ytemp(lats)
62  rlatmin = MIN( ytemp(1), ytemp(lats) )
63  rlatmax = MAX( ytemp(1), ytemp(lats) )
64     
65  IF     (rlatmin>=-pis2-0.5.AND.rlatmax<=pis2+0.5) THEN; radlat = .TRUE.
66  ELSE IF(rlatmin>= -90.-0.5.AND.rlatmax<= 90.+0.5) THEN; radlat = .FALSE.
67  ELSE; WRITE(6,*) 'Problems with data latitudes for file '//TRIM(title)
68  END IF
69
70  IF(.NOT.radlon) xtemp(:)=xtemp(:)*pi/180.
71  IF(.NOT.radlat) ytemp(:)=ytemp(:)*pi/180.
72
73!--- FLIPPED LONGITUDES
74  IF(invlon) THEN
75    champf(:,:)=champd(:,:); xf(:)=xtemp(:)
76
77  !--- Longitudes rotated to get them in [ -pi pi] interval.
78    DO i=1,lons;      IF(xf(i)>pi) EXIT;                   END DO; ip180 = i
79    DO i=1,lons;      IF(xf(i)>pi) xf(i)=xf(i)-depi;       END DO
80    DO i= ip180,lons; ind=i-ip180+1; xtemp(ind)=xf(i);     END DO
81    DO i= ind+1,lons;                xtemp(i  )=xf(i-ind); END DO
82
83  !--- Longitudes rotated in champf
84    DO j=1,lats
85      DO i=ip180,lons; ind=i-ip180+1; champd(ind,j)=champf(i    ,j); END DO
86      DO i=ind+1,lons;                champd(i  ,j)=champf(i-ind,j); END DO
87    END DO
88  END IF
89
90!--- FLIPPED LATITUDES
91  IF(invlat) THEN
92    yf(:)=ytemp(:)
93    champf(:,:)=champd(:,:)
94    ytemp(lats:1:-1)=yf(:)
95    DO j=1,lats; champd(:,lats-j+1)=champf(:,j); END DO
96  END IF
97
98!--- FOR BARYCENTRIC INTERPOLATION
99  IF(interbar) THEN
100    oldxd1 = xtemp(1)
101    xtemp(1:lons-1)=0.5*(xtemp(1:lons-1)+xtemp(2:lons))
102    xtemp(  lons  )=0.5*(xtemp(  lons  )+oldxd1+depi)
103    ytemp(1:lats-1)=0.5*(ytemp(1:lats-1)+ytemp(2:lats))
104  END IF
105  DEALLOCATE(champf)
106  xf(:)=xtemp(:); DEALLOCATE(xtemp)
107  yf(:)=ytemp(:); DEALLOCATE(ytemp)
108
109END SUBROUTINE conf_dat2d
110!
111!-------------------------------------------------------------------------------
112
113
114!-------------------------------------------------------------------------------
115!
116SUBROUTINE conf_dat3d(title, xd, yd, zd, xf, yf, zf, champd, interbar)
117!
118!-------------------------------------------------------------------------------
119! Author: P. Le Van
120!-------------------------------------------------------------------------------
121! Purpose: Configure the 3D data field "champd" so that:
122!     - Longitudes are in [ -pi     pi   ]
123!     - Latitudes  are in [  pi/2. -pi/2 ]
124!     - Vertical levels from ground to model top (in Pascals)
125!  * xd / yd are initial lon / lats.
126!  * xf / yf are output  lon / lats, possibly modified to satisfy configuration.
127!  * zf      are output pressures
128!  * interbar is TRUE for barycentric interpolation.
129!-------------------------------------------------------------------------------
130  USE assert_eq_m, ONLY: assert_eq
131  IMPLICIT NONE
132!-------------------------------------------------------------------------------
133! Arguments:
134  CHARACTER(LEN=*), INTENT(IN)    :: title
135  REAL,             INTENT(IN)    :: xd(:), yd(:), zd(:) ! (lons) (lats) (levs)
136  REAL,             INTENT(INOUT) :: xf(:), yf(:), zf(:) ! (lons) (lats) (levs)
137  REAL,             INTENT(INOUT) :: champd(:,:,:)       ! (lons,lats,levs)
138  LOGICAL,          INTENT(IN)    :: interbar
139!-------------------------------------------------------------------------------
140! Local variables:
141  CHARACTER(LEN=256) :: modname="conf_dat3d"
142  INTEGER :: i, j, l, ip180, ind, lons, lats, levs
143  LOGICAL :: radlon, invlon ,radlat, invlat, invlev
144  REAL    :: pi, pis2, depi, presmax, rlatmin, rlatmax, oldxd1
145  REAL, ALLOCATABLE :: xtemp(:) , ytemp(:), ztemp(:), champf(:,:,:)
146!-------------------------------------------------------------------------------
147  lons=assert_eq(SIZE(xd),SIZE(xf),SIZE(champd,1),TRIM(modname)//" lons")
148  lats=assert_eq(SIZE(yd),SIZE(yf),SIZE(champd,2),TRIM(modname)//" lats")
149  levs=assert_eq(SIZE(zd),SIZE(zf),SIZE(champd,3),TRIM(modname)//" levs")
150  ALLOCATE(xtemp(lons),ytemp(lats),ztemp(levs),champf(lons,lats,levs))
151  xtemp(:)=xd(:); ytemp(:)=yd(:); ztemp(:)=zd(:)
152  pi = 2. * ASIN(1.)
153  pis2 = pi/2.
154  depi = 2. * pi
155  radlon=.FALSE.; invlon=.FALSE.
156  IF     (xtemp(1)>=-pi-0.5.AND.xtemp(lons)<=  pi+0.5) THEN
157    radlon = .TRUE.;  invlon = .FALSE.
158  ELSE IF(xtemp(1)>=   -0.5.AND.xtemp(lons)<=depi+0.5) THEN
159    radlon = .TRUE.;  invlon = .TRUE.
160  ELSE IF(xtemp(1)>= -180.5.AND.xtemp(lons)<=   180.5) THEN
161    radlon = .FALSE.; invlon = .FALSE.
162  ELSE IF(xtemp(1)>=   -0.5.AND.xtemp(lons)<=   360.5) THEN
163    radlon = .FALSE.; invlon = .TRUE.
164  ELSE; WRITE(6,*) 'Problems with data longitudes for file '//TRIM(title)
165  END IF
166  invlat = ytemp(1)<ytemp(lats)
167  rlatmin = MIN( ytemp(1), ytemp(lats) )
168  rlatmax = MAX( ytemp(1), ytemp(lats) )
169     
170  IF     (rlatmin>=-pis2-0.5.AND.rlatmax<=pis2+0.5) THEN; radlat = .TRUE.
171  ELSE IF(rlatmin>= -90.-0.5.AND.rlatmax<= 90.+0.5) THEN; radlat = .FALSE.
172  ELSE; WRITE(6,*) 'Problems with data latitudes for file '//TRIM(title)
173  END IF
174
175  IF(.NOT.radlon) xtemp(:)=xtemp(:)*pi/180.
176  IF(.NOT.radlat) ytemp(:)=ytemp(:)*pi/180.
177
178!--- FLIPPED LONGITUDES
179  IF(invlon) THEN
180    champf(:,:,:)=champd(:,:,:); xf(:)=xtemp(:)
181
182  !--- Longitudes rotated to get them in [ -pi pi] interval.
183    DO i=1,lons;      IF(xf(i)>pi) EXIT;                   END DO; ip180 = i
184    DO i=1,lons;      IF(xf(i)>pi) xf(i)=xf(i)-depi;       END DO
185    DO i= ip180,lons; ind=i-ip180+1; xtemp(ind)=xf(i);     END DO
186    DO i= ind+1,lons;                xtemp(i  )=xf(i-ind); END DO
187
188  !--- Longitudes rotated in champf
189    DO l=1,levs
190    DO j=1,lats
191      DO i=ip180,lons; ind=i-ip180+1; champd(ind,j,l)=champf(i    ,j,l); END DO
192      DO i=ind+1,lons;                champd(i  ,j,l)=champf(i-ind,j,l); END DO
193    END DO
194    END DO
195  END IF
196
197!--- FLIPPED LATITUDES
198  IF(invlat) THEN
199    yf(:)=ytemp(:)
200    champf(:,:,:)=champd(:,:,:)
201    ytemp(lats:1:-1)=yf(:)
202    DO l=1,levs
203      DO j=1,lats; champd(:,lats-j+1,l)=champf(:,j,l); END DO
204    END DO
205  END IF
206
207!--- FOR BARYCENTRIC INTERPOLATION
208  IF(interbar) THEN
209    oldxd1 = xtemp(1)
210    xtemp(1:lons-1)=0.5*(xtemp(1:lons-1)+xtemp(2:lons))
211    xtemp(  lons  )=0.5*(xtemp(  lons  )+oldxd1+depi)
212    ytemp(1:lats-1)=0.5*(ytemp(1:lats-1)+ytemp(2:lats))
213  END IF
214
215!--- FLIPPED LEVELS
216  invlev=ztemp(1)<ztemp(levs)
217  IF(MAX(ztemp(1),ztemp(levs))<1200.) ztemp(:)=ztemp(:)*100.
218  IF(invlev) THEN
219    zf(:)=ztemp(:)
220    champf(:,:,:)=champd(:,:,:)
221    ztemp(levs:1:-1)=zf(:)
222    DO l=1,levs; champd(:,:,levs+1-l)=champf(:,:,l); END DO
223  END IF
224
225  DEALLOCATE(champf)
226  xf(:)=xtemp(:); DEALLOCATE(xtemp)
227  yf(:)=ytemp(:); DEALLOCATE(ytemp)
228  zf(:)=ztemp(:); DEALLOCATE(ztemp)
229
230END SUBROUTINE conf_dat3d
231!
232!-------------------------------------------------------------------------------
233
234END MODULE conf_dat_m
235
Note: See TracBrowser for help on using the repository browser.