
 Fujitsu Fortran Version 3.0  Tue Oct  8 20:52:53 2002

 Compilation information
   Current directory : /home/takepiro/work/spmodel/spml-0.1/src
   Source file       : wghi_module.f90

 Option information
   Command line options : -c -I/usr/lib/gtool4-ifc6/include -I/usr/lib/dcl-5.2-ifc6/include -O -Qoption,ld,-Bdynamic -Vaxlib -w -Am

 Module "wghi_module"
  (line-no.)(nest)
          1        !----------------------------------------------------------------------
          2        !     Copyright (c) 2002 Shin-ichi Takehiro. All rights reserved.
          3        !----------------------------------------------------------------------
          4        !ɽ  wghi_module
          5        !
          6        !  2002/05/26  ݹ
          7        !
          8        !ľؤΤȤ
          9        !
         10        !   -------------     nz         g(nz)
         11        !		     nz-1/2      h(nz)
         12        !   -------------     nz-1       g(nz-1), i(nz-1)
         13        !   ...
         14        !
         15        !   ...
         16        !   -------------     2          g(2), i(2)
         17        !		     3/2         h(2)
         18        !   -------------     1          g(1), i(1)
         19        !		     1/2         h(1)
         20        !   -------------  iz=0          g(0)
         21        !
         22        !
         23        !      ǡ index
         24        !
         25        !   x   : ٳʻ
         26        !   y   : ٳʻ
         27        !   g   : ľʻ()      Grid points
         28        !   h   : ľʻ(Ⱦ)    grid points with Half integer
         29        !   i   : ľʻ(,ΰ) Internal grid points
         30        !
         31        !   w   : ʿĴȡ
         32        !   a   : Ǥդμ
         33        !
         34        !   xyg, xyh : 3 ʻǡ
         35        !   wg, wh   : ʿڥȥľʻǡ(ľʻ)
         36        !
         37        !   xy      : ʿʻǡ
         38        !   w       : ʿڥȥǡ
         39        !
         40        !   xg, xh  : ٱ̳ʻǡ
         41        !   yg, yh  : Ҹ̳ʻǡ
         42        !
         43        module wghi_module
         44        
         45          use wa_module
         46          implicit none
         47          private
         48        
         49          public wghi_initial
         50        
         51         !---- ɸѿ ----
         52          public x_lon, x_lon_weight
         53          public y_lat, y_lat_weight
         54          public g_ver, g_ver_weight
         55          public h_ver, h_ver_weight
         56          public i_ver, i_ver_weight
         57          public l_nm, nm_l
         58          public xyg_lon, xyg_lat, xyg_ver
         59          public xyh_lon, xyh_lat, xyh_ver
         60          public xyi_lon, xyi_lat, xyi_ver
         61        
         62         !---- ǡѴ ----
         63          public w_xy, xy_w                                     ! ڥȥѴ
         64          public xyg_wg, xyh_wh, xyi_wi                         ! ڥȥѴ
         65          public wg_xyg, wh_xyh, wi_xyi                         ! ڥȥѴ
         66          public xyi_xyg, xi_xg, yi_yg, i_g, wi_wg              ! ľʻ
         67          public xyi_xyh, xi_xh, yi_yh, i_h, wi_wh              ! ľʻ
         68          public xyh_xyg, xh_xg, yh_yg, h_g, wh_wg              ! ľʻ
         69          public xyg_yg, xyh_yh, xyi_yi                         ! Ҹ̥ǡĥ
         70        
         71         !---- ʬ ----
         72        
         73          public wi_dver_wh, wh_dver_wg                            ! ľʬ
         74          public xyg_gradlon_wg, xyg_gradlat_wg                    ! ʿ
         75          public xyh_gradlon_wh, xyh_gradlat_wh                    ! ʿ
         76          public xyi_gradlon_wi, xyi_gradlat_wi                    ! ʿ
         77          public wg_divlon_xyg, wg_divlat_xyg                      ! ȯ
         78          public wh_divlon_xyh, wh_divlat_xyh, wh_div_xyh_xyh_xyg  ! ȯ
         79          public wi_divlon_xyi, wi_divlat_xyi, wi_div_xyi_xyi_xyh  ! ȯ
         80          public wg_laplah_wg, wh_laplah_wh, wi_laplah_wi          ! ʿץ饷
         81        
         82         !---- ʬ ----
         83          public yg_intlon_xyg, yh_intlon_xyh
         84          public xg_intlat_xyg, xh_intlat_xyh
         85          public xy_intver_xyg, xy_intver_xyh
         86          public y_intlonver_xyg, y_intlonver_xyh
         87          public x_intlatver_xyg, x_intlatver_xyh
         88          public g_intlonlat_xyg, h_intlonlat_xyh
         89          public intlonlatver_xyg, intlonlatver_xyh
         90        
         91         !---- ʿ ----
         92          public yg_avrlon_xyg, yh_avrlon_xyh
         93          public xg_avrlat_xyg, xh_avrlat_xyh
         94          public xy_avrver_xyg, xy_avrver_xyh
         95          public x_avrlatver_xyg, x_avrlatver_xyh
         96          public y_avrlonver_xyg, y_avrlonver_xyh
         97          public g_avrlonlat_xyg, h_avrlonlat_xyh
         98          public avrlonlatver_xyg, avrlonlatver_xyh
         99        
        100        
        101         !---- ѿ ----
        102          integer            :: im=64, jm=32, km=16  ! ʻ(, , ľ)
        103          integer            :: nm=21                ! ʿȿ
        104          real(8)            :: ra=1.0               ! Ⱦ
        105          real(8), parameter :: pi=3.1415926535897932385D0
        106        
        107          real(8), dimension(:), allocatable :: g_ver, g_ver_weight    ! ľɸ
        108          real(8), dimension(:), allocatable :: h_ver, h_ver_weight    ! ľɸ
        109          real(8), dimension(:), allocatable :: i_ver, i_ver_weight    ! ľɸ
        110        
        111          real(8), dimension(:,:,:), allocatable :: xyg_lon, xyg_lat, xyg_ver ! ɸ
        112          real(8), dimension(:,:,:), allocatable :: xyh_lon, xyh_lat, xyh_ver ! ɸ
        113          real(8), dimension(:,:,:), allocatable :: xyi_lon, xyi_lat, xyi_ver ! ɸ
        114        
        115          save im, jm, km, nm, ra
        116        
        117          contains
        118          !---------------  -----------------
        119            subroutine wghi_initial(i,j,k,n,r,g_in,h_in)
        120        
        121             integer,intent(in) :: i, j, k        ! ʻ(, , ľ)
        122             integer,intent(in) :: n              ! ȿ(ʿ)
        123             real(8),intent(in) :: r              ! Ⱦ
        124             real(8),intent(in),dimension(0:k) :: g_in       ! ľɸ
        125             real(8),intent(in),dimension(k)   :: h_in       ! ľɸ
        126        
        127             im = i  ; jm = j ; km = k
        128             nm = n  ; ra = r
        129        
        130             call wa_initial(nm,im,jm,km+1)
        131        
        132             allocate(g_ver(0:km))
        133             allocate(g_ver_weight(0:km))
        134             allocate(h_ver(km))
        135             allocate(h_ver_weight(km))
        136             allocate(i_ver(km-1))
        137             allocate(i_ver_weight(km-1))
        138        
        139             g_ver = g_in ; h_ver = h_in
        140        
        141             h_ver_weight = g_ver(1:km) - g_ver(0:km-1)
        142             g_ver_weight(1:km-1) = h_ver(2:km) - h_ver(1:km-1)
        143             g_ver_weight(0) = h_ver(1) - g_ver(0)
        144             g_ver_weight(km) = g_ver(km) - h_ver(km)
        145        
        146             i_ver = i_g(g_ver) ; i_ver_weight = i_g(g_ver_weight)
        147        
        148             allocate(xyg_lon(im,jm,0:km))
        149             allocate(xyg_lat(im,jm,0:km))
        150             allocate(xyg_ver(im,jm,0:km))
        151             allocate(xyh_lon(im,jm,km))
        152             allocate(xyh_lat(im,jm,km))
        153             allocate(xyh_ver(im,jm,km))
        154             allocate(xyi_lon(im,jm,km-1))
        155             allocate(xyi_lat(im,jm,km-1))
        156             allocate(xyi_ver(im,jm,km-1))
        157        
        158             xyg_lon = spread(xy_lon,3,km+1)
        159             xyg_lat = spread(xy_lat,3,km+1)
        160             xyg_ver = spread(spread(g_ver,1,jm),1,im)
        161             xyh_lon = spread(xy_lon,3,km)
        162             xyh_lat = spread(xy_lat,3,km)
        163             xyh_ver = spread(spread(h_ver,1,jm),1,im)
        164             xyi_lon = spread(xy_lon,3,km-1)
        165             xyi_lat = spread(xy_lat,3,km-1)
        166             xyi_ver = spread(spread(i_ver,1,jm),1,im)
        167        
        168           end subroutine wghi_initial
        169        
        170          !--------------- Ѵ -----------------
        171        
        172            function xyg_wg(wg_data)  ! ڥȥ -> ʻ
        173              real(8), dimension(im,jm,0:km)                     :: xyg_wg
        174              real(8), dimension((nm+1)*(nm+1),0:km), intent(in) :: wg_data
        175        
        176              xyg_wg = xya_wa(wg_data)
        177        
        178            end function xyg_wg
        179        
        180            function xyh_wh(wh_data)  ! ڥȥ -> ʻ
        181              real(8), dimension(im,jm,km)                     :: xyh_wh
        182              real(8), dimension((nm+1)*(nm+1),km), intent(in) :: wh_data
        183        
        184              xyh_wh = xya_wa(wh_data)
        185        
        186            end function xyh_wh
        187        
        188            function xyi_wi(wi_data)  ! ڥȥ -> ʻ
        189              real(8), dimension(im,jm,km-1)                     :: xyi_wi
        190              real(8), dimension((nm+1)*(nm+1),km-1), intent(in) :: wi_data
        191        
        192              xyi_wi = xya_wa(wi_data)
        193        
        194            end function xyi_wi
        195        
        196            function wg_xyg(xyg_data)  ! ʻ -> ڥȥ
        197              real(8), dimension((nm+1)*(nm+1),0:km)             :: wg_xyg
        198              real(8), dimension(im,jm,0:km), intent(in)         :: xyg_data
        199        
        200              wg_xyg = wa_xya(xyg_data)
        201        
        202            end function wg_xyg
        203        
        204            function wh_xyh(xyh_data)  ! ʻ -> ڥȥ
        205              real(8), dimension((nm+1)*(nm+1),km)             :: wh_xyh
        206              real(8), dimension(im,jm,km), intent(in)         :: xyh_data
        207        
        208              wh_xyh = wa_xya(xyh_data)
        209        
        210            end function wh_xyh
        211        
        212            function wi_xyi(xyi_data)  ! ʻ -> ڥȥ
        213              real(8), dimension((nm+1)*(nm+1),km-1)             :: wi_xyi
        214              real(8), dimension(im,jm,km-1), intent(in)         :: xyi_data
        215        
        216              wi_xyi = wa_xya(xyi_data)
        217        
        218            end function wi_xyi
        219        
        220            function xyi_xyg(xyg_data) ! ľѴ
        221              real(8), dimension(im,jm,km-1)                     :: xyi_xyg
        222              real(8), dimension(im,jm,0:km), intent(in)         :: xyg_data
        223        
        224              xyi_xyg = xyg_data(:,:,1:km-1)
        225            end function xyi_xyg
        226        
        227            function xi_xg(xg_data)    ! ľѴ
        228              real(8), dimension(im,km-1)                     :: xi_xg
        229              real(8), dimension(im,0:km), intent(in)         :: xg_data
        230        
        231              xi_xg = xg_data(:,1:km-1)
        232            end function xi_xg
        233        
        234            function yi_yg(yg_data)    ! ľѴ
        235              real(8), dimension(jm,km-1)                     :: yi_yg
        236              real(8), dimension(jm,0:km), intent(in)         :: yg_data
        237        
        238              yi_yg = yg_data(:,1:km-1)
        239            end function yi_yg
        240        
        241            function i_g(g_data)       ! ľѴ
        242              real(8), dimension(km-1)                     :: i_g
        243              real(8), dimension(0:km), intent(in)         :: g_data
        244        
        245              i_g = g_data(1:km-1)
        246            end function i_g
        247        
        248            function wi_wg(wg_data)    ! ľѴ
        249              real(8), dimension((nm+1)*(nm+1),km-1)             :: wi_wg
        250              real(8), dimension((nm+1)*(nm+1),0:km), intent(in) :: wg_data
        251        
        252              wi_wg = wg_data(:,1:km-1)
        253            end function wi_wg
        254        
        255            function xyi_xyh(xyh_data)    ! ľ(ľ)
        256              real(8), dimension(im,jm,km-1)                   :: xyi_xyh
        257              real(8), dimension(im,jm,km), intent(in)         :: xyh_data
        258              integer :: k
        259        
        260     1        do k=1,km-1
        261     1           xyi_xyh(:,:,k) = (   (h_ver(k+1)-i_ver(k))*xyh_data(:,:,k)  &
        262     1                              + (i_ver(k)-h_ver(k))*xyh_data(:,:,k+1) )&
        263     1                            /(h_ver(k+1)-h_ver(k))
        264     1        enddo
        265            end function xyi_xyh
        266        
        267            function xi_xh(xh_data)    ! ľ(ľ)
        268              real(8), dimension(im,km-1)                   :: xi_xh
        269              real(8), dimension(im,km), intent(in)         :: xh_data
        270              integer :: k
        271        
        272     1        do k=1,km-1
        273     1           xi_xh(:,k) = (   (h_ver(k+1)-i_ver(k))*xh_data(:,k)  &
        274     1                          + (i_ver(k)-h_ver(k))*xh_data(:,k+1) )&
        275     1                            /(h_ver(k+1)-h_ver(k))
        276     1        enddo
        277            end function xi_xh
        278        
        279            function yi_yh(yh_data)    ! ľ(ľ)
        280              real(8), dimension(jm,km-1)                   :: yi_yh
        281              real(8), dimension(jm,km), intent(in)         :: yh_data
        282              integer :: k
        283        
        284     1        do k=1,km-1
        285     1           yi_yh(:,k) = (   (h_ver(k+1)-i_ver(k))*yh_data(:,k)  &
        286     1                          + (i_ver(k)-h_ver(k))*yh_data(:,k+1) )&
        287     1                            /(h_ver(k+1)-h_ver(k))
        288     1        enddo
        289            end function yi_yh
        290        
        291            function i_h(h_data)       ! ľ(ľ)
        292              real(8), dimension(km-1)                   :: i_h
        293              real(8), dimension(km), intent(in)         :: h_data
        294              integer :: k
        295        
        296     1        do k=1,km-1
        297     1           i_h(k) = (   (h_ver(k+1)-i_ver(k))*h_data(k)  &
        298     1                      + (i_ver(k)-h_ver(k))*h_data(k+1) )&
        299     1                    /(h_ver(k+1)-h_ver(k))
        300     1        enddo
        301            end function i_h
        302        
        303            function wi_wh(wh_data)    ! ľ(ľ)
        304              real(8), dimension((nm+1)*(nm+1),km-1)           :: wi_wh
        305              real(8), dimension((nm+1)*(nm+1),km), intent(in) :: wh_data
        306              integer :: k
        307        
        308     1        do k=1,km-1
        309     1           wi_wh(:,k) = (   (h_ver(k+1)-i_ver(k))*wh_data(:,k)  &
        310     1                          + (i_ver(k)-h_ver(k))*wh_data(:,k+1) )&
        311     1                            /(h_ver(k+1)-h_ver(k))
        312     1        enddo
        313            end function wi_wh
        314        
        315            function xyh_xyg(xyg_data)    ! ľ(ľ)
        316              real(8), dimension(im,jm,km)                   :: xyh_xyg
        317              real(8), dimension(im,jm,0:km), intent(in)     :: xyg_data
        318              integer :: k
        319        
        320     1        do k=1,km
        321     1           xyh_xyg(:,:,k) = (   (g_ver(k)-h_ver(k))*xyg_data(:,:,k-1)  &
        322     1                              + (h_ver(k)-g_ver(k-1))*xyg_data(:,:,k) )&
        323     1                            /(g_ver(k)-g_ver(k-1))
        324     1        enddo
        325            end function xyh_xyg
        326        
        327            function xh_xg(xg_data)    ! ľ(ľ)
        328              real(8), dimension(im,km)                       :: xh_xg
        329              real(8), dimension(im,0:km), intent(in)         :: xg_data
        330              integer :: k
        331        
        332     1        do k=1,km
        333     1           xh_xg(:,k) = (   (g_ver(k)-h_ver(k))*xg_data(:,k-1)  &
        334     1                          + (h_ver(k)-g_ver(k-1))*xg_data(:,k) )&
        335     1                          /(g_ver(k)-g_ver(k-1))
        336     1        enddo
        337            end function xh_xg
        338        
        339            function yh_yg(yg_data)    ! ľ(ľ)
        340              real(8), dimension(jm,km)                       :: yh_yg
        341              real(8), dimension(jm,0:km), intent(in)         :: yg_data
        342              integer :: k
        343        
        344     1        do k=1,km
        345     1           yh_yg(:,k) = (   (g_ver(k)-h_ver(k))*yg_data(:,k-1)  &
        346     1                          + (h_ver(k)-g_ver(k-1))*yg_data(:,k) )&
        347     1                          /(g_ver(k)-g_ver(k-1))
        348     1        enddo
        349            end function yh_yg
        350        
        351            function h_g(g_data)       ! ľ(ľ)
        352              real(8), dimension(km)                       :: h_g
        353              real(8), dimension(0:km), intent(in)         :: g_data
        354              integer :: k
        355        
        356     1        do k=1,km
        357     1           h_g(k) = (   (g_ver(k)-h_ver(k))*g_data(k-1)  &
        358     1                      + (h_ver(k)-g_ver(k-1))*g_data(k) ) &
        359     1                     /(g_ver(k)-g_ver(k-1))
        360     1        enddo
        361            end function h_g
        362        
        363            function wh_wg(wg_data)    ! ľ(ľ)
        364              real(8), dimension((nm+1)*(nm+1),km)               :: wh_wg
        365              real(8), dimension((nm+1)*(nm+1),0:km), intent(in) :: wg_data
        366              integer :: k
        367        
        368     1        do k=1,km
        369     1           wh_wg(:,k) = (   (g_ver(k)-h_ver(k))*wg_data(:,k-1)  &
        370     1                          + (h_ver(k)-g_ver(k-1))*wg_data(:,k) )&
        371     1                          /(g_ver(k)-g_ver(k-1))
        372     1        enddo
        373            end function wh_wg
        374        
        375            function xyg_yg(yg_data)    ! Ҹ̥ǡĥ
        376              real(8), dimension(im,jm,0:km)              :: xyg_yg
        377              real(8), dimension(jm,0:km), intent(in)     :: yg_data
        378        
        379              xyg_yg = spread(yg_data,1,im)
        380            end function xyg_yg
        381        
        382            function xyh_yh(yh_data)    ! Ҹ̥ǡĥ
        383              real(8), dimension(im,jm,km)                :: xyh_yh
        384              real(8), dimension(jm,km), intent(in)       :: yh_data
        385        
        386              xyh_yh = spread(yh_data,1,im)
        387            end function xyh_yh
        388        
        389            function xyi_yi(yi_data)    ! Ҹ̥ǡĥ
        390              real(8), dimension(im,jm,km-1)              :: xyi_yi
        391              real(8), dimension(jm,km-1), intent(in)     :: yi_data
        392        
        393              xyi_yi = spread(yi_data,1,im)
        394            end function xyi_yi
        395        
        396          !--------------- ʬ׻ -----------------
        397            function wi_dver_wh(wh_data)            ! ľʬ
        398              real(8), dimension((nm+1)*(nm+1),km-1), intent(in) :: wh_data
        399              real(8), dimension((nm+1)*(nm+1),0:km)             :: wi_dver_wh
        400        
        401              wi_dver_wh = (wh_data(:,2:km) - wh_data(:,1:km-1)) &
        402                                      /spread(i_ver_weight,1,(nm+1)*(nm+1))
        403        
        404            end function wi_dver_wh
        405        
        406            function wh_dver_wg(wg_data)            ! ľʬ
        407              real(8), dimension((nm+1)*(nm+1),0:km), intent(in) :: wg_data
        408              real(8), dimension((nm+1)*(nm+1),km)               :: wh_dver_wg
        409        
        410              wh_dver_wg = (wg_data(:,1:km) - wg_data(:,0:km-1)) &
        411                            /spread(h_ver_weight,1,(nm+1)*(nm+1))
        412        
        413            end function wh_dver_wg
        414        
        415            function xyg_gradlon_wg(wg_data) ! ڥȥ˺Ѥ۷ʬ
        416                                             ! 1/rcosա/ߦ
        417              real(8), dimension((nm+1)*(nm+1),0:km), intent(in) :: wg_data
        418              real(8), dimension(im,jm,0:km)                     :: xyg_gradlon_wg
        419        
        420              xyg_gradlon_wg = xya_gradlon_wa(wg_data)/ra
        421        
        422            end function xyg_gradlon_wg
        423        
        424            function xyh_gradlon_wh(wh_data) ! ڥȥ˺Ѥ۷ʬ
        425                                             ! 1/rcosա/ߦ
        426              real(8), dimension((nm+1)*(nm+1),km), intent(in) :: wh_data
        427              real(8), dimension(im,jm,km)                     :: xyh_gradlon_wh
        428        
        429              xyh_gradlon_wh = xya_gradlon_wa(wh_data)/ra
        430        
        431            end function xyh_gradlon_wh
        432        
        433            function xyi_gradlon_wi(wi_data) ! ڥȥ˺Ѥ۷ʬ
        434                                             ! 1/rcosա/ߦ
        435              real(8), dimension((nm+1)*(nm+1),km-1), intent(in) :: wi_data
        436              real(8), dimension(im,jm,km-1)                     :: xyi_gradlon_wi
        437        
        438              xyi_gradlon_wi = xya_gradlon_wa(wi_data)/ra
        439        
        440            end function xyi_gradlon_wi
        441        
        442            function xyg_gradlat_wg(wg_data) ! ڥȥ˺Ѥ۷ʬ
        443                                             ! 1/r /ߦ
        444              real(8), dimension((nm+1)*(nm+1),0:km), intent(in) :: wg_data
        445              real(8), dimension(im,jm,0:km)                     :: xyg_gradlat_wg
        446        
        447              xyg_gradlat_wg = xya_gradlat_wa(wg_data)/ra
        448            end function xyg_gradlat_wg
        449        
        450            function xyh_gradlat_wh(wh_data) ! ڥȥ˺Ѥ۷ʬ
        451                                             ! 1/r /ߦ
        452              real(8), dimension((nm+1)*(nm+1),km), intent(in) :: wh_data
        453              real(8), dimension(im,jm,km)                     :: xyh_gradlat_wh
        454        
        455              xyh_gradlat_wh = xya_gradlat_wa(wh_data)/ra
        456            end function xyh_gradlat_wh
        457        
        458            function xyi_gradlat_wi(wi_data) ! ڥȥ˺Ѥ۷ʬ
        459                                             ! 1/r /ߦ
        460              real(8), dimension((nm+1)*(nm+1),km-1), intent(in) :: wi_data
        461              real(8), dimension(im,jm,km-1)                     :: xyi_gradlat_wi
        462        
        463              xyi_gradlat_wi = xya_gradlat_wa(wi_data)/ra
        464            end function xyi_gradlat_wi
        465        
        466            function wg_divlon_xyg(xyg_data)   ! ʻҤ˺Ѥȯʬ
        467                                               ! 1/rcosա/ߦ
        468        
        469              real(8), dimension(im,jm,0:km), intent(in)   :: xyg_data
        470              real(8), dimension((nm+1)*(nm+1),0:km)       :: wg_divlon_xyg
        471        
        472              wg_divlon_xyg = wa_divlon_xya(xyg_data)/ra
        473            end function wg_divlon_xyg
        474        
        475            function wh_divlon_xyh(xyh_data)   ! ʻҤ˺Ѥȯʬ
        476                                               ! 1/rcosա/ߦ
        477        
        478              real(8), dimension(im,jm,km), intent(in)   :: xyh_data
        479              real(8), dimension((nm+1)*(nm+1),km)       :: wh_divlon_xyh
        480        
        481              wh_divlon_xyh = wa_divlon_xya(xyh_data)/ra
        482            end function wh_divlon_xyh
        483        
        484            function wi_divlon_xyi(xyi_data)   ! ʻҤ˺Ѥȯʬ
        485                                               ! 1/rcosա/ߦ
        486        
        487              real(8), dimension(im,jm,km-1), intent(in)   :: xyi_data
        488              real(8), dimension((nm+1)*(nm+1),km-1)       :: wi_divlon_xyi
        489        
        490              wi_divlon_xyi = wa_divlon_xya(xyi_data)/ra
        491            end function wi_divlon_xyi
        492        
        493            function wg_divlat_xyg(xyg_data)   ! ʻҤ˺Ѥȯʬ
        494                                               ! 1/rcosա(f cos)/ߦ
        495        
        496              real(8), dimension(im,jm,0:km), intent(in)   :: xyg_data
        497              real(8), dimension((nm+1)*(nm+1),0:km)       :: wg_divlat_xyg
        498        
        499              wg_divlat_xyg = wa_divlat_xya(xyg_data)/ra
        500            end function wg_divlat_xyg
        501        
        502            function wh_divlat_xyh(xyh_data)   ! ʻҤ˺Ѥȯʬ
        503                                               ! 1/rcosա(f cos)/ߦ
        504        
        505              real(8), dimension(im,jm,km), intent(in)   :: xyh_data
        506              real(8), dimension((nm+1)*(nm+1),km)       :: wh_divlat_xyh
        507        
        508              wh_divlat_xyh = wa_divlat_xya(xyh_data)/ra
        509            end function wh_divlat_xyh
        510        
        511            function wi_divlat_xyi(xyi_data)   ! ʻҤ˺Ѥȯʬ
        512                                               ! 1/rcosա(f cos)/ߦ
        513        
        514              real(8), dimension(im,jm,km-1), intent(in)   :: xyi_data
        515              real(8), dimension((nm+1)*(nm+1),km-1)       :: wi_divlat_xyi
        516        
        517              wi_divlat_xyi = wa_divlat_xya(xyi_data)/ra
        518            end function wi_divlat_xyi
        519        
        520            function wi_div_xyi_xyi_xyh(xyi_vlon,xyi_vlat,xyh_vver)    ! ȯ
        521        
        522              real(8), dimension(im,jm,km-1), intent(in) :: xyi_vlon   ! ʬ
        523              real(8), dimension(im,jm,km-1), intent(in) :: xyi_vlat   ! ʬ
        524              real(8), dimension(im,jm,km), intent(in)   :: xyh_vver   ! ưʬ
        525              real(8), dimension((nm+1)*(nm+1),km-1)     :: wi_div_xyi_xyi_xyh
        526        
        527              wi_div_xyi_xyi_xyh =   wi_divlon_xyi(xyi_vlon) &
        528                                   + wi_divlat_xyi(xyi_vlat) &
        529                                   + wi_dver_wh(wh_xyh(xyh_vver))
        530        
        531            end function wi_div_xyi_xyi_xyh
        532        
        533            function wh_div_xyh_xyh_xyg(xyh_vlon,xyh_vlat,xyg_vver)    ! ȯ
        534        
        535              real(8), dimension(im,jm,km), intent(in)   :: xyh_vlon   ! ʬ
        536              real(8), dimension(im,jm,km), intent(in)   :: xyh_vlat   ! ʬ
        537              real(8), dimension(im,jm,0:km), intent(in) :: xyg_vver   ! ưʬ
        538              real(8), dimension((nm+1)*(nm+1),km)       :: wh_div_xyh_xyh_xyg
        539        
        540              wh_div_xyh_xyh_xyg =   wh_divlon_xyh(xyh_vlon) &
        541                                   + wh_divlat_xyh(xyh_vlat) &
        542                                   + wh_dver_wg(wg_xyg(xyg_vver))
        543        
        544            end function wh_div_xyh_xyh_xyg
        545        
        546            function wg_laplah_wg(wg_data) ! ڥȥ˺Ѥʿץ饷
        547                                           ! ^2_H
        548              real(8), dimension((nm+1)*(nm+1),0:km), intent(in) :: wg_data
        549              real(8), dimension((nm+1)*(nm+1),0:km)             :: wg_laplah_wg
        550        
        551              wg_laplah_wg = wa_lapla_wa(wg_data)
        552        
        553            end function wg_laplah_wg
        554        
        555            function wh_laplah_wh(wh_data) ! ڥȥ˺Ѥʿץ饷
        556                                           ! ^2_H
        557              real(8), dimension((nm+1)*(nm+1),km), intent(in) :: wh_data
        558              real(8), dimension((nm+1)*(nm+1),km)             :: wh_laplah_wh
        559        
        560              wh_laplah_wh = wa_lapla_wa(wh_data)
        561        
        562            end function wh_laplah_wh
        563        
        564            function wi_laplah_wi(wi_data) ! ڥȥ˺Ѥʿץ饷
        565                                           ! ^2_H
        566              real(8), dimension((nm+1)*(nm+1),km-1), intent(in) :: wi_data
        567              real(8), dimension((nm+1)*(nm+1),km-1)             :: wi_laplah_wi
        568        
        569              wi_laplah_wi = wa_lapla_wa(wi_data)
        570        
        571            end function wi_laplah_wi
        572        
        573          !--------------- ʬ׻ -----------------
        574        
        575            function yg_intlon_xyg(xyg_data)  ! (Ӿ)ʬ
        576              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data ! 3 ʻ
        577              real(8), dimension(jm,0:km)  :: yg_intlon_xyg          ! Ҹ̳ʻ
        578        
        579              yg_intlon_xyg = ya_intlon_xya(xyg_data)
        580            end function yg_intlon_xyg
        581        
        582            function yh_intlon_xyh(xyh_data)  ! (Ӿ)ʬ
        583              real(8), dimension(im,jm,km), intent(in) :: xyh_data ! 3 ʻ
        584              real(8), dimension(jm,km)  :: yh_intlon_xyh          ! Ҹ̳ʻ
        585        
        586              yh_intlon_xyh = ya_intlon_xya(xyh_data)
        587            end function yh_intlon_xyh
        588        
        589            function xg_intlat_xyg(xyg_data)  ! ʬ
        590              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data ! 3 ʻ
        591              real(8), dimension(im,0:km)  :: xg_intlat_xyg          ! ٱ߳ʻ
        592        
        593              xg_intlat_xyg = xa_intlat_xya(xyg_data)
        594            end function xg_intlat_xyg
        595        
        596            function xh_intlat_xyh(xyh_data)  ! ʬ
        597              real(8), dimension(im,jm,km), intent(in) :: xyh_data ! 3 ʻ
        598              real(8), dimension(im,km)  :: xh_intlat_xyh          ! ٱ߳ʻ
        599        
        600              xh_intlat_xyh = xa_intlat_xya(xyh_data)
        601            end function xh_intlat_xyh
        602        
        603            function xy_intver_xyg(xyg_data)  ! ưʬ
        604              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data ! 3 ʻ
        605              real(8), dimension(im,jm)  :: xy_intver_xyg            ! ʿʻ
        606              integer :: k
        607        
        608              xy_intver_xyg = 0
        609     1        do k=0,km
        610     1           xy_intver_xyg(:,:) = xy_intver_xyg(:,:) &
        611     1                         + xyg_data(:,:,k) * g_ver_weight(k)
        612     1        enddo
        613            end function xy_intver_xyg
        614        
        615            function xy_intver_xyh(xyh_data)  ! ưʬ
        616              real(8), dimension(im,jm,km), intent(in) :: xyh_data ! 3 ʻ
        617              real(8), dimension(im,jm)  :: xy_intver_xyh            ! ʿʻ
        618              integer :: k
        619        
        620              xy_intver_xyh = 0
        621     1        do k=1,km
        622     1           xy_intver_xyh(:,:) = xy_intver_xyh(:,:) &
        623     1                         + xyh_data(:,:,k) * h_ver_weight(k)
        624     1        enddo
        625            end function xy_intver_xyh
        626        
        627            function x_intlatver_xyg(xyg_data)  ! ư(Ҹ)ʬ
        628              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data  ! 3 ʻ
        629              real(8), dimension(im)     :: x_intlatver_xyg      ! ٳʻ
        630        
        631              x_intlatver_xyg = x_intlat_xy(xy_intver_xyg(xyg_data))
        632            end function x_intlatver_xyg
        633        
        634            function x_intlatver_xyh(xyh_data)  ! ư(Ҹ)ʬ
        635              real(8), dimension(im,jm,km), intent(in) :: xyh_data  ! 3 ʻ
        636              real(8), dimension(im)     :: x_intlatver_xyh         ! ٳʻ
        637        
        638              x_intlatver_xyh = x_intlat_xy(xy_intver_xyh(xyh_data))
        639            end function x_intlatver_xyh
        640        
        641            function y_intlonver_xyg(xyg_data)  ! ư(ٱ)ʬ
        642              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data  ! 3 ʻ
        643              real(8), dimension(jm)       :: y_intlonver_xyg         ! ٳʻ
        644        
        645              y_intlonver_xyg = y_intlon_xy(xy_intver_xyg(xyg_data))
        646            end function y_intlonver_xyg
        647        
        648            function y_intlonver_xyh(xyh_data)  ! ư(ٱ)ʬ
        649              real(8), dimension(im,jm,km), intent(in) :: xyh_data  ! 3 ʻ
        650              real(8), dimension(jm)       :: y_intlonver_xyh       ! ٳʻ
        651        
        652              y_intlonver_xyh = y_intlon_xy(xy_intver_xyh(xyh_data))
        653            end function y_intlonver_xyh
        654        
        655            function g_intlonlat_xyg(xyg_data)  ! ٷ(ʿ)ʬ
        656              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data  ! 3 ʻ
        657              real(8), dimension(0:km)     :: g_intlonlat_xyg         ! ư³ʻ
        658        
        659              g_intlonlat_xyg = a_intlonlat_xya(xyg_data)
        660            end function g_intlonlat_xyg
        661        
        662            function h_intlonlat_xyh(xyh_data)  ! ٷ(ʿ)ʬ
        663              real(8), dimension(im,jm,km), intent(in) :: xyh_data  ! 3 ʻ
        664              real(8), dimension(km)     :: h_intlonlat_xyh         ! ư³ʻ
        665        
        666              h_intlonlat_xyh = a_intlonlat_xya(xyh_data)
        667            end function h_intlonlat_xyh
        668        
        669            function intlonlatver_xyg(xyg_data) ! ٷư()ʬ
        670              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data  ! 3 ʻ
        671              real(8)                        :: intlonlatver_xyg      ! ʬ
        672        
        673              intlonlatver_xyg = intlonlat_xy(xy_intver_xyg(xyg_data))
        674            end function intlonlatver_xyg
        675        
        676            function intlonlatver_xyh(xyh_data) ! ٷư()ʬ
        677              real(8), dimension(im,jm,km), intent(in) :: xyh_data  ! 3 ʻ
        678              real(8)                        :: intlonlatver_xyh      ! ʬ
        679        
        680              intlonlatver_xyh = intlonlat_xy(xy_intver_xyh(xyh_data))
        681            end function intlonlatver_xyh
        682        
        683          !--------------- ʿѷ׻ -----------------
        684        
        685            function yg_avrlon_xyg(xyg_data)  ! (Ӿ)ʬ
        686              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data ! 3 ʻ
        687              real(8), dimension(jm,0:km)  :: yg_avrlon_xyg          ! Ҹ̳ʻ
        688        
        689              yg_avrlon_xyg = ya_avrlon_xya(xyg_data)
        690            end function yg_avrlon_xyg
        691        
        692            function yh_avrlon_xyh(xyh_data)  ! (Ӿ)ʬ
        693              real(8), dimension(im,jm,km), intent(in) :: xyh_data ! 3 ʻ
        694              real(8), dimension(jm,km)  :: yh_avrlon_xyh          ! Ҹ̳ʻ
        695        
        696              yh_avrlon_xyh = ya_avrlon_xya(xyh_data)
        697            end function yh_avrlon_xyh
        698        
        699            function xg_avrlat_xyg(xyg_data)  ! ʬ
        700              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data ! 3 ʻ
        701              real(8), dimension(im,0:km)  :: xg_avrlat_xyg          ! ٱ߳ʻ
        702        
        703              xg_avrlat_xyg = xa_avrlat_xya(xyg_data)
        704            end function xg_avrlat_xyg
        705        
        706            function xh_avrlat_xyh(xyh_data)  ! ʬ
        707              real(8), dimension(im,jm,km), intent(in) :: xyh_data ! 3 ʻ
        708              real(8), dimension(im,km)  :: xh_avrlat_xyh          ! ٱ߳ʻ
        709        
        710              xh_avrlat_xyh = xa_avrlat_xya(xyh_data)
        711            end function xh_avrlat_xyh
        712        
        713            function xy_avrver_xyg(xyg_data)  ! ưʬ
        714              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data ! 3 ʻ
        715              real(8), dimension(im,jm)  :: xy_avrver_xyg            ! ʿʻ
        716        
        717              xy_avrver_xyg = xy_intver_xyg(xyg_data)/sum(g_ver_weight)
        718            end function xy_avrver_xyg
        719        
        720            function xy_avrver_xyh(xyh_data)  ! ưʬ
        721              real(8), dimension(im,jm,km), intent(in) :: xyh_data ! 3 ʻ
        722              real(8), dimension(im,jm)  :: xy_avrver_xyh            ! ʿʻ
        723        
        724              xy_avrver_xyh = xy_intver_xyh(xyh_data)/sum(h_ver_weight)
        725            end function xy_avrver_xyh
        726        
        727            function x_avrlatver_xyg(xyg_data)  ! ư(Ҹ)ʬ
        728              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data  ! 3 ʻ
        729              real(8), dimension(im)     :: x_avrlatver_xyg      ! ٳʻ
        730        
        731              x_avrlatver_xyg = x_avrlat_xy(xy_avrver_xyg(xyg_data))
        732            end function x_avrlatver_xyg
        733        
        734            function x_avrlatver_xyh(xyh_data)  ! ư(Ҹ)ʬ
        735              real(8), dimension(im,jm,km), intent(in) :: xyh_data  ! 3 ʻ
        736              real(8), dimension(im)     :: x_avrlatver_xyh         ! ٳʻ
        737        
        738              x_avrlatver_xyh = x_avrlat_xy(xy_avrver_xyh(xyh_data))
        739            end function x_avrlatver_xyh
        740        
        741            function y_avrlonver_xyg(xyg_data)  ! ư(ٱ)ʬ
        742              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data  ! 3 ʻ
        743              real(8), dimension(jm)       :: y_avrlonver_xyg         ! ٳʻ
        744        
        745              y_avrlonver_xyg = y_avrlon_xy(xy_avrver_xyg(xyg_data))
        746            end function y_avrlonver_xyg
        747        
        748            function y_avrlonver_xyh(xyh_data)  ! ư(ٱ)ʬ
        749              real(8), dimension(im,jm,km), intent(in) :: xyh_data  ! 3 ʻ
        750              real(8), dimension(jm)       :: y_avrlonver_xyh       ! ٳʻ
        751        
        752              y_avrlonver_xyh = y_avrlon_xy(xy_avrver_xyh(xyh_data))
        753            end function y_avrlonver_xyh
        754        
        755            function g_avrlonlat_xyg(xyg_data)  ! ٷ(ʿ)ʬ
        756              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data  ! 3 ʻ
        757              real(8), dimension(0:km)     :: g_avrlonlat_xyg         ! ư³ʻ
        758        
        759              g_avrlonlat_xyg = a_avrlonlat_xya(xyg_data)
        760            end function g_avrlonlat_xyg
        761        
        762            function h_avrlonlat_xyh(xyh_data)  ! ٷ(ʿ)ʬ
        763              real(8), dimension(im,jm,km), intent(in) :: xyh_data  ! 3 ʻ
        764              real(8), dimension(km)     :: h_avrlonlat_xyh         ! ư³ʻ
        765        
        766              h_avrlonlat_xyh = a_avrlonlat_xya(xyh_data)
        767            end function h_avrlonlat_xyh
        768        
        769            function avrlonlatver_xyg(xyg_data) ! ٷư()ʬ
        770              real(8), dimension(im,jm,0:km), intent(in) :: xyg_data  ! 3 ʻ
        771              real(8)                        :: avrlonlatver_xyg      ! ʬ
        772        
        773              avrlonlatver_xyg = avrlonlat_xy(xy_avrver_xyg(xyg_data))
        774            end function avrlonlatver_xyg
        775        
        776            function avrlonlatver_xyh(xyh_data) ! ٷư()ʬ
        777              real(8), dimension(im,jm,km), intent(in) :: xyh_data  ! 3 ʻ
        778              real(8)                        :: avrlonlatver_xyh      ! ʬ
        779        
        780              avrlonlatver_xyh = avrlonlat_xy(xy_avrver_xyh(xyh_data))
        781            end function avrlonlatver_xyh
        782        
        783          end module wghi_module

 Procedure information
   Lines      : 783
   Statements : 498

 Total information
   Procedures       : 1
   Total lines      : 783
   Total statements : 498

