/[lmdze]/trunk/filtrez/inifilr.f
ViewVC logotype

Diff of /trunk/filtrez/inifilr.f

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

revision 168 by guez, Wed Sep 9 10:41:47 2015 UTC revision 169 by guez, Mon Sep 14 17:13:16 2015 UTC
# Line 2  module inifilr_m Line 2  module inifilr_m
2    
3    IMPLICIT NONE    IMPLICIT NONE
4    
   ! North:  
   
5    INTEGER jfiltnu, jfiltnv    INTEGER jfiltnu, jfiltnv
6    ! index of the last scalar line filtered in northern hemisphere    ! index of the last line filtered in northern hemisphere at rlat[uv]
7      ! latitudes
   real, pointer:: matriceun(:, :, :) ! (iim, iim, jfiltnu - 1)  
   ! matrice filtre pour les champs situes sur la grille scalaire  
   
   real, pointer:: matrinvn(:, :, :) ! (iim, iim, jfiltnu - 1)  
   ! matrice filtre pour les champs situes sur la grille scalaire, pour  
   ! le filtre inverse  
   
   real, pointer:: matricevn(:, :, :) ! (iim, iim, jfiltnv)  
   ! matrice filtre pour les champs situes sur la grille de V ou de Z  
   
   ! South:  
8    
9    integer jfiltsu, jfiltsv    integer jfiltsu, jfiltsv
10    ! index of the first line filtered in southern hemisphere    ! index of the first line filtered in southern hemisphere at
11      ! rlat[uv] latitudes
12    
13      ! Filtre pour les champs situes sur la grille scalaire (longitudes
14      ! rlonv, latitudes rlatu) :
15      real, pointer:: matriceun(:, :, :) ! (iim, iim, jfiltnu - 1)
16    real, pointer:: matriceus(:, :, :) ! (iim, iim, jjm - jfiltsu + 1)    real, pointer:: matriceus(:, :, :) ! (iim, iim, jjm - jfiltsu + 1)
   ! matrice filtre pour les champs situes sur la grille scalaire  
17    
18      ! Filtre pour les champs situes sur la grille scalaire (longitudes
19      ! rlonv, latitudes rlatu), pour le filtre inverse :
20      real, pointer:: matrinvn(:, :, :) ! (iim, iim, jfiltnu - 1)
21    real, pointer:: matrinvs(:, :, :) ! (iim, iim, jjm - jfiltsu + 1)    real, pointer:: matrinvs(:, :, :) ! (iim, iim, jjm - jfiltsu + 1)
   ! matrice filtre pour les champs situes sur la grille scalaire, pour  
   ! le filtre inverse  
22    
23      ! Filtre pour les champs situes sur la grille de la vorticit\'e
24      ! (longitudes rlonu, latitudes rlatv)
25      real, pointer:: matricevn(:, :, :) ! (iim, iim, jfiltnv) matrice
26    real, pointer:: matricevs(:, :, :) ! (iim, iim, jjm - jfiltsv + 1)    real, pointer:: matricevs(:, :, :) ! (iim, iim, jjm - jfiltsv + 1)
   ! matrice filtre pour les champs situes sur la grille de V ou de Z  
27    
28  contains  contains
29    
# Line 51  contains Line 44  contains
44      use inifgn_m, only: inifgn      use inifgn_m, only: inifgn
45      use inifilr_hemisph_m, only: inifilr_hemisph      use inifilr_hemisph_m, only: inifilr_hemisph
46      use jumble, only: new_unit      use jumble, only: new_unit
47      use nr_util, only: pi, ifirstloc      use nr_util, only: pi, ifirstloc, assert
48    
49      ! Local:      ! Local:
50    
# Line 60  contains Line 53  contains
53      real eignvl(iim) ! eigenvalues (<= 0) sorted in descending order      real eignvl(iim) ! eigenvalues (<= 0) sorted in descending order
54      INTEGER j, unit      INTEGER j, unit
55      REAL colat0 ! > 0      REAL colat0 ! > 0
56      integer j1 ! index of smallest positive latitude      integer j1 ! index of negative latitude closest to the equator
57    
58      real eignfnu(iim, iim), eignfnv(iim, iim)      real eignfnu(iim, iim), eignfnv(iim, iim)
59      ! eigenvectors of the discrete second derivative with respect to longitude      ! eigenvectors of the discrete second derivative with respect to
60        ! longitude, at rlon[uv] longitudes
61    
62      !-----------------------------------------------------------      !-----------------------------------------------------------
63    
# Line 76  contains Line 70  contains
70      colat0 = min(0.5, minval(dlatu) / minval(xprimu(:iim)))      colat0 = min(0.5, minval(dlatu) / minval(xprimu(:iim)))
71      PRINT *, 'colat0 = ', colat0      PRINT *, 'colat0 = ', colat0
72    
73      rlamda = iim / (pi * colat0 / grossismx) / sqrt(- eignvl(2: iim))      rlamda = iim / pi / colat0 * grossismx / sqrt(- eignvl(2: iim))
74      print *, "1 / rlamda(iim) = ", 1. / rlamda(iim)      print *, "1 / rlamda(iim) = ", 1. / rlamda(iim)
75    
76        ! This is demonstrated in the notes but just to be sure:
77        call assert(rlamda(iim) * colat0 >= 1. - epsilon(0.), &
78             "inifilr rlamda(iim) * colat0")
79    
80      call new_unit(unit)      call new_unit(unit)
81      open(unit, file = "modfrst.csv", status = "replace", action = "write")      open(unit, file = "modfrst.csv", status = "replace", action = "write")
82      write(unit, fmt = *) '"rlat (degrees)" modfrst' ! title line      write(unit, fmt = *) '"rlat (degrees)" modfrst' ! title line
83    
84      j1 = ifirstloc(rlatu <= 0.)      j1 = ifirstloc(rlatu <= 0.)
85    
86      call inifilr_hemisph(rlatu(j1 - 1:2:- 1), colat0, rlamda, unit, eignfnv, &      call inifilr_hemisph(rlatu(j1 - 1:2:- 1), rlamda, unit, eignfnv, jfiltnu, &
87           jfiltnu, matriceun, matrinvn)           matriceun, matrinvn)
88      jfiltnu = j1 - jfiltnu      jfiltnu = j1 - jfiltnu
89      matriceun = matriceun(:, :, jfiltnu - 1:1:- 1)      matriceun = matriceun(:, :, jfiltnu - 1:1:- 1)
90      matrinvn = matrinvn(:, :, jfiltnu - 1:1:- 1)      matrinvn = matrinvn(:, :, jfiltnu - 1:1:- 1)
91    
92      call inifilr_hemisph(- rlatu(j1:jjm), colat0, rlamda, unit, eignfnv, &      call inifilr_hemisph(- rlatu(j1:jjm), rlamda, unit, eignfnv, jfiltsu, &
93           jfiltsu, matriceus, matrinvs)           matriceus, matrinvs)
94      jfiltsu = j1 - 1 + jfiltsu      jfiltsu = j1 - 1 + jfiltsu
95    
96      j1 = ifirstloc(rlatv <= 0.)      j1 = ifirstloc(rlatv <= 0.)
97    
98      call inifilr_hemisph(rlatv(j1 - 1:1:- 1), colat0, rlamda, unit, eignfnu, &      call inifilr_hemisph(rlatv(j1 - 1:1:- 1), rlamda, unit, eignfnu, jfiltnv, &
99           jfiltnv, matricevn)           matricevn)
100      jfiltnv = j1 - jfiltnv      jfiltnv = j1 - jfiltnv
101      matricevn = matricevn(:, :, jfiltnv:1:- 1)      matricevn = matricevn(:, :, jfiltnv:1:- 1)
102    
103      call inifilr_hemisph(- rlatv(j1:jjm), colat0, rlamda, unit, eignfnu, &      call inifilr_hemisph(- rlatv(j1:jjm), rlamda, unit, eignfnu, jfiltsv, &
104           jfiltsv, matricevs)           matricevs)
105      jfiltsv = j1 - 1 + jfiltsv      jfiltsv = j1 - 1 + jfiltsv
106    
107      close(unit)      close(unit)
Legend:
Removed from v.168  
changed lines
  Added in v.169
  ViewVC Help
Powered by ViewVC 1.1.21