/[lmdze]/trunk/phylmd/Orography/grid_noro_m.f
ViewVC logotype

Diff of /trunk/phylmd/Orography/grid_noro_m.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 134 by guez, Wed Apr 29 15:47:56 2015 UTC revision 147 by guez, Wed Jun 17 14:20:14 2015 UTC
# Line 9  contains Line 9  contains
9    
10      ! From dyn3d/grid_noro.F, version 1.1.1.1 2004/05/19 12:53:06      ! From dyn3d/grid_noro.F, version 1.1.1.1 2004/05/19 12:53:06
11    
12      ! Authors: François Lott, Laurent Li, A. Harzallah and Laurent      ! Authors: Fran\c{}cois Lott, Laurent Li, A. Harzallah and Laurent
13      ! Fairhead      ! Fairhead
14    
15      ! Compute the parameters of the sub-grid scale orography scheme as      ! Compute the parameters of the sub-grid scale orography scheme as
# Line 32  contains Line 32  contains
32      use mva9_m, only: mva9      use mva9_m, only: mva9
33      use nr_util, only: assert, pi      use nr_util, only: assert, pi
34    
35      REAL, intent(in):: xdata(:), ydata(:) ! coordinates of input field      ! Coordinates of input field:
36      REAL, intent(in):: zdata(:, :) ! input field      REAL, intent(in):: xdata(:) ! (iusn)
37        REAL, intent(in):: ydata(:) ! (jusn)
38    
39        REAL, intent(in):: zdata(:, :) ! (iusn, jusn) input field
40      REAL, intent(in):: x(:), y(:) ! coordinates of output field      REAL, intent(in):: x(:), y(:) ! coordinates of output field
41    
42      ! Correlations of US Navy orography gradients:      ! Correlations of US Navy orography gradients:
43      REAL, intent(out):: zphi(:, :) ! orography not smoothed      REAL, intent(out):: zphi(:, :) ! (iim + 1, jjm + 1) orography not smoothed
44      real, intent(out):: zmea(:, :) ! smoothed orography      real, intent(out):: zmea(:, :) ! (iim + 1, jjm + 1) smoothed orography
45      real, intent(out):: zstd(:, :) ! Standard deviation      real, intent(out):: zstd(:, :) ! (iim + 1, jjm + 1) Standard deviation
46      REAL, intent(out):: zsig(:, :) ! Slope      REAL, intent(out):: zsig(:, :) ! (iim + 1, jjm + 1) Slope
47      real, intent(out):: zgam(:, :) ! Anisotropy      real, intent(out):: zgam(:, :) ! (iim + 1, jjm + 1) Anisotropy
48      real, intent(out):: zthe(:, :) ! Orientation of the small axis  
49      REAL, intent(out):: zpic(:, :) ! Maximum altitude      real, intent(out):: zthe(:, :) ! (iim + 1, jjm + 1)
50      real, intent(out):: zval(:, :) ! Minimum altitude      ! Orientation of the small axis
51    
52      real, intent(out):: mask(:, :) ! fraction of land      REAL, intent(out):: zpic(:, :) ! (iim + 1, jjm + 1) Maximum altitude
53        real, intent(out):: zval(:, :) ! (iim + 1, jjm + 1) Minimum altitude
54    
55        real, intent(out):: mask(:, :) ! (iim + 1, jjm + 1) fraction of land
56    
57      ! Variables local to the procedure:      ! Variables local to the procedure:
58    
# Line 63  contains Line 69  contains
69      ! Correlations of US Navy orography gradients:      ! Correlations of US Navy orography gradients:
70      REAL, dimension(iusn + 2 * iext, jusn + 2):: zxtzxusn, zytzyusn, zxtzyusn      REAL, dimension(iusn + 2 * iext, jusn + 2):: zxtzxusn, zytzyusn, zxtzyusn
71    
72      real, dimension(iim + 1, jjm + 1):: mask_tmp, num_tot, num_lan, zmea0      real, dimension(iim + 1, jjm + 1):: mask_tmp, num_tot, num_lan
73      REAL a(iim + 1), b(iim + 1), c(jjm + 1), d(jjm + 1)      REAL a(iim + 1), b(iim + 1), c(jjm + 1), d(jjm + 1)
74      real weighx, weighy, xincr, xk, xp, xm, xw, xq, xl      real weighx, weighy, xincr, xk, xp, xm, xw, xq, xl
75      real zbordnor, zdeltax, zbordsud, zdeltay, zbordoue, zlenx, zleny, zmeasud      real zbordnor, zdeltax, zbordsud, zdeltay, zbordoue, zlenx, zleny, zmeasud
76      real zllmpic, zllmmea, zllmgam, zllmthe, zllmstd, zllmsig, zllmval      real zweinor, zweisud, zmeanor, zbordest
     real zpicnor, zminthe, zsigsud, zstdnor, zstdsud, zvalsud, zvalnor  
     real zweinor, zweisud, zsignor, zpicsud, zmeanor, zbordest  
77      integer ii, i, jj, j      integer ii, i, jj, j
78      real, parameter:: rad = 6371229.      real, parameter:: rad = 6371229.
79    
# Line 88  contains Line 92  contains
92           size(zsig, 2), size(zgam, 2), size(zthe, 2), size(zpic, 2), &           size(zsig, 2), size(zgam, 2), size(zthe, 2), size(zpic, 2), &
93           size(zval, 2), size(mask, 2)/) == jjm + 1, "grid_noro jjm")           size(zval, 2), size(mask, 2)/) == jjm + 1, "grid_noro jjm")
94    
95      print *, "Paramètres de l'orographie à l'échelle sous-maille"      print *, "Parameters of subgrid-scale orography"
96      zdeltay = 2. * pi / real(jusn) * rad      zdeltay = 2. * pi / real(jusn) * rad
97    
98      ! Extension of the US Navy database for computations at boundaries:      ! Extension of the US Navy database for computations at boundaries:
# Line 216  contains Line 220  contains
220      ! Compute parameters needed by the Lott & Miller (1997) and Lott      ! Compute parameters needed by the Lott & Miller (1997) and Lott
221      ! (1999) subgrid-scale orographic scheme.      ! (1999) subgrid-scale orographic scheme.
222    
     zllmmea = 0.  
     zllmstd = 0.  
     zllmsig = 0.  
     zllmgam = 0.  
     zllmpic = 0.  
     zllmval = 0.  
     zllmthe = 0.  
     zminthe = 0.  
223      DO ii = 1, iim + 1      DO ii = 1, iim + 1
224         DO jj = 1, jjm + 1         DO jj = 1, jjm + 1
225            mask(ii, jj) = num_lan(ii, jj) / num_tot(ii, jj)            mask(ii, jj) = num_lan(ii, jj) / num_tot(ii, jj)
# Line 248  contains Line 244  contains
244         zytzy(ii, jjm + 1) = zytzy(ii, jjm)         zytzy(ii, jjm + 1) = zytzy(ii, jjm)
245      ENDDO      ENDDO
246    
247      zmea0 = zmea ! not smoothed      ! Masque prenant en compte maximum de terre. On met un seuil \`a 10
248        ! % de terre car en dessous les param\`etres de surface n'ont pas de
249        ! sens.
250        mask_tmp = merge(1., 0., mask >= 0.1)
251    
252        zphi(:iim, :) = zmea(:iim, :) * mask_tmp(:iim, :)
253        ! (zmea is not yet smoothed)
254    
255      ! Filters to smooth out fields for input into subgrid-scale      ! Filters to smooth out fields for input into subgrid-scale
256      ! orographic scheme.      ! orographic scheme.
# Line 262  contains Line 264  contains
264      CALL MVA9(zxtzy)      CALL MVA9(zxtzy)
265      CALL MVA9(zytzy)      CALL MVA9(zytzy)
266    
     ! Masque prenant en compte maximum de terre. On met un seuil à 10  
     ! % de terre car en dessous les paramètres de surface n'ont pas de  
     ! sens.  
     mask_tmp = merge(1., 0., mask >= 0.1)  
   
267      DO ii = 1, iim      DO ii = 1, iim
268         DO jj = 1, jjm + 1         DO jj = 1, jjm + 1
269            ! Coefficients K, L et M:            ! Coefficients K, L et M:
# Line 276  contains Line 273  contains
273            xp = xk - sqrt(xl**2 + xm**2)            xp = xk - sqrt(xl**2 + xm**2)
274            xq = xk + sqrt(xl**2 + xm**2)            xq = xk + sqrt(xl**2 + xm**2)
275            xw = 1e-8            xw = 1e-8
276            if(xp.le.xw) xp = 0.            if (xp <= xw) xp = 0.
277            if(xq.le.xw) xq = xw            if (xq <= xw) xq = xw
278            if(abs(xm).le.xw) xm = xw * sign(1., xm)            if (abs(xm) <= xw) xm = xw * sign(1., xm)
279            ! modification pour masque de terre fractionnaire            ! modification pour masque de terre fractionnaire
280            ! slope:            ! slope:
281            zsig(ii, jj) = sqrt(xq) * mask_tmp(ii, jj)            zsig(ii, jj) = sqrt(xq) * mask_tmp(ii, jj)
# Line 286  contains Line 283  contains
283            zgam(ii, jj) = xp / xq * mask_tmp(ii, jj)            zgam(ii, jj) = xp / xq * mask_tmp(ii, jj)
284            ! angle theta:            ! angle theta:
285            zthe(ii, jj) = 57.29577951 * atan2(xm, xl) / 2. * mask_tmp(ii, jj)            zthe(ii, jj) = 57.29577951 * atan2(xm, xl) / 2. * mask_tmp(ii, jj)
           zphi(ii, jj) = zmea0(ii, jj) * mask_tmp(ii, jj)  
           zmea(ii, jj) = zmea(ii, jj) * mask_tmp(ii, jj)  
           zpic(ii, jj) = zpic(ii, jj) * mask_tmp(ii, jj)  
           zval(ii, jj) = zval(ii, jj) * mask_tmp(ii, jj)  
           zstd(ii, jj) = zstd(ii, jj) * mask_tmp(ii, jj)  
           zllmmea = MAX(zmea(ii, jj), zllmmea)  
           zllmstd = MAX(zstd(ii, jj), zllmstd)  
           zllmsig = MAX(zsig(ii, jj), zllmsig)  
           zllmgam = MAX(zgam(ii, jj), zllmgam)  
           zllmthe = MAX(zthe(ii, jj), zllmthe)  
           zminthe = min(zthe(ii, jj), zminthe)  
           zllmpic = MAX(zpic(ii, jj), zllmpic)  
           zllmval = MAX(zval(ii, jj), zllmval)  
286         ENDDO         ENDDO
287      ENDDO      ENDDO
288    
289      print *, 'MEAN ORO: ', zllmmea      zmea(:iim, :) = zmea(:iim, :) * mask_tmp(:iim, :)
290      print *, 'ST. DEV.: ', zllmstd      zpic(:iim, :) = zpic(:iim, :) * mask_tmp(:iim, :)
291      print *, 'PENTE: ', zllmsig      zval(:iim, :) = zval(:iim, :) * mask_tmp(:iim, :)
292      print *, 'ANISOTROP: ', zllmgam      zstd(:iim, :) = zstd(:iim, :) * mask_tmp(:iim, :)
293      print *, 'ANGLE: ', zminthe, zllmthe  
294      print *, 'pic: ', zllmpic      print *, 'MEAN ORO: ', MAXVAL(zmea(:iim, :))
295      print *, 'val: ', zllmval      print *, 'ST. DEV.: ', MAXVAL(zstd(:iim, :))
296        print *, 'PENTE: ', MAXVAL(zsig(:iim, :))
297        print *, 'ANISOTROP: ', MAXVAL(zgam(:iim, :))
298        print *, 'ANGLE: ', minval(zthe(:iim, :)), MAXVAL(zthe(:iim, :))
299        print *, 'pic: ', MAXVAL(zpic(:iim, :))
300        print *, 'val: ', MAXVAL(zval(:iim, :))
301    
302      ! gamma and theta at 1. and 0. at poles      ! gamma and theta at 1. and 0. at poles
303      zmea(iim + 1, :) = zmea(1, :)      zmea(iim + 1, :) = zmea(1, :)
# Line 320  contains Line 309  contains
309      zgam(iim + 1, :) = zgam(1, :)      zgam(iim + 1, :) = zgam(1, :)
310      zthe(iim + 1, :) = zthe(1, :)      zthe(iim + 1, :) = zthe(1, :)
311    
312      zmeanor = 0.      zweinor = sum(weight(:iim, 1))
313      zmeasud = 0.      zweisud = sum(weight(:iim, jjm + 1))
314      zstdnor = 0.      zmeanor = sum(zmea(:iim, 1) * weight(:iim, 1))
315      zstdsud = 0.      zmeasud = sum(zmea(:iim, jjm + 1) * weight(:iim, jjm + 1))
     zsignor = 0.  
     zsigsud = 0.  
     zweinor = 0.  
     zweisud = 0.  
     zpicnor = 0.  
     zpicsud = 0.  
     zvalnor = 0.  
     zvalsud = 0.  
   
     DO ii = 1, iim  
        zweinor = zweinor + weight(ii, 1)  
        zweisud = zweisud + weight(ii, jjm + 1)  
        zmeanor = zmeanor + zmea(ii, 1) * weight(ii, 1)  
        zmeasud = zmeasud + zmea(ii, jjm + 1) * weight(ii, jjm + 1)  
        zstdnor = zstdnor + zstd(ii, 1) * weight(ii, 1)  
        zstdsud = zstdsud + zstd(ii, jjm + 1) * weight(ii, jjm + 1)  
        zsignor = zsignor + zsig(ii, 1) * weight(ii, 1)  
        zsigsud = zsigsud + zsig(ii, jjm + 1) * weight(ii, jjm + 1)  
        zpicnor = zpicnor + zpic(ii, 1) * weight(ii, 1)  
        zpicsud = zpicsud + zpic(ii, jjm + 1) * weight(ii, jjm + 1)  
        zvalnor = zvalnor + zval(ii, 1) * weight(ii, 1)  
        zvalsud = zvalsud + zval(ii, jjm + 1) * weight(ii, jjm + 1)  
     ENDDO  
316    
317      zmea(:, 1) = zmeanor / zweinor      zmea(:, 1) = zmeanor / zweinor
318      zmea(:, jjm + 1) = zmeasud / zweisud      zmea(:, jjm + 1) = zmeasud / zweisud
# Line 354  contains Line 320  contains
320      zphi(:, 1) = zmeanor / zweinor      zphi(:, 1) = zmeanor / zweinor
321      zphi(:, jjm + 1) = zmeasud / zweisud      zphi(:, jjm + 1) = zmeasud / zweisud
322    
323      zpic(:, 1) = zpicnor / zweinor      zpic(:, 1) = sum(zpic(:iim, 1) * weight(:iim, 1)) / zweinor
324      zpic(:, jjm + 1) = zpicsud / zweisud      zpic(:, jjm + 1) = sum(zpic(:iim, jjm + 1) * weight(:iim, jjm + 1)) &
325             / zweisud
326      zval(:, 1) = zvalnor / zweinor  
327      zval(:, jjm + 1) = zvalsud / zweisud      zval(:, 1) = sum(zval(:iim, 1) * weight(:iim, 1)) / zweinor
328        zval(:, jjm + 1) = sum(zval(:iim, jjm + 1) * weight(:iim, jjm + 1)) &
329      zstd(:, 1) = zstdnor / zweinor           / zweisud
330      zstd(:, jjm + 1) = zstdsud / zweisud  
331        zstd(:, 1) = sum(zstd(:iim, 1) * weight(:iim, 1)) / zweinor
332      zsig(:, 1) = zsignor / zweinor      zstd(:, jjm + 1) = sum(zstd(:iim, jjm + 1) * weight(:iim, jjm + 1)) &
333      zsig(:, jjm + 1) = zsigsud / zweisud           / zweisud
334    
335        zsig(:, 1) = sum(zsig(:iim, 1) * weight(:iim, 1)) / zweinor
336        zsig(:, jjm + 1) = sum(zsig(:iim, jjm + 1) * weight(:iim, jjm + 1)) &
337             / zweisud
338    
339      zgam(:, 1) = 1.      zgam(:, 1) = 1.
340      zgam(:, jjm + 1) = 1.      zgam(:, jjm + 1) = 1.
Legend:
Removed from v.134  
changed lines
  Added in v.147
  ViewVC Help
Powered by ViewVC 1.1.21