!= dcrtm 驩帷峩繝輔Λ繝�け繧ｹ 險育ｮ� 繝励Ο繧ｰ繝ｩ繝�
!
!= dcrtm vertical flux calculation main program
!
! Authors::   Masanori Onishi
! Version::   $Id: main_dcrtm.f90,v 1.00 onishi$
! Tag Name::  $Name: dcrtm-20160705 $
! Copyright:: Copyright (C) GFD Dennou Club, 2008-2015. All rights reserved.
! License::   See COPYRIGHT[link:../../../COPYRIGHT]

program sum_verticalflux
  !
  ! <b>Note that Japanese and English are described in parallel.</b>
  !

  ! 繝｢繧ｸ繝･繝ｼ繝ｫ蠑慕畑 ; USE statements
  !

  ! netcdf 繝��繧ｿ蜈･蜃ｺ蜉�
  ! netcdf data input or output
  !
  use gtool_history

  ! 螳｣險譁� ; Declaration statements
  !
  implicit none

  ! 險ｭ螳壼､画焚, 莉�
  ! setteing parameter, values
  !
  character(99) :: InFileName
                          ! 蜈･蜉� 繝輔ぃ繧､繝ｫ蜷�
                          ! input file name
  character(99) :: OutFileName
                          ! 蜃ｺ蜉帙ヵ繧｡繧､繝ｫ蜷�
                          ! output file name
  real(8) :: ResWN
                          ! 豕｢謨ｰ蛻�ｧ｣閭ｽ [m-1]
                          ! resolution of wavenumber [m-1]
  real(8), allocatable :: r_Press(:)
                          ! 豌怜悸 (蜊頑紛謨ｰ繝ｬ繝吶Ν) [Pa]
                          ! pressure (half level) [Pa]
  real(8), allocatable :: rw_RadFluxUp(:,:)
                          ! 諠第弌螟ｧ豌励�蟆��縺ｫ逕ｱ譚･縺吶ｋ荳雁髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2 m]
                          ! upward flux originated from emittance of planetary atmospheire [W m-2 m]
  real(8), allocatable :: rw_RadFluxDn(:,:)
                          ! 諠第弌螟ｧ豌励�蟆��縺ｫ逕ｱ譚･縺吶ｋ荳句髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2 m]
                          ! downward flux originated from emittance of planetary atmospheire [W m-2 m]
  real(8), allocatable :: rw_SolFluxUp(:,:)
                          ! 荳ｭ蠢�弌蜈･蟆�↓逕ｱ譚･縺吶ｋ荳雁髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2 m]
                          ! upward flux originated from emittance of central star [W m-2 m]
  real(8), allocatable :: rw_SolFluxDn(:,:)
                          ! 荳ｭ蠢�弌蜈･蟆�↓逕ｱ譚･縺吶ｋ荳句髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2 m]
                          ! downward flux originated from emittance of central star [W m-2 m]

  integer :: fnml = 12
                          ! 陬�ｽｮ逡ｪ蜿ｷ(name list file)
                          ! device number : namelist file
  integer :: wmax
                          ! 豕｢謨ｰ繝ｬ繝吶Ν譛螟ｧ蛟､
                          ! wavenumber level maximum
  integer :: kmax
                          ! 驩帷峩繝ｬ繝吶Ν譛螟ｧ蛟､
                          ! vertical level maximum
  integer :: k_w, k_r
                          ! 菴懈･ｭ螟画焚
                          ! work variables
  real(8), pointer :: a_Axis(:) => null()
                          ! 驟榊�縺ｮ繧ｵ繧､繧ｺ繧定ｪｿ縺ｹ繧九◆繧√�驟榊�
                          ! dummy array

  ! 險育ｮ励＆繧後ｋ迚ｩ逅�㍼
  ! calculated values
  !
  real(8), allocatable :: r_RadFluxUp(:)
                          ! 諠第弌螟ｧ豌励�蟆��縺ｫ逕ｱ譚･縺吶ｋ荳雁髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2]
                          ! upward flux originated from emittance of planetary atmospheire [W m-2]
  real(8), allocatable :: r_RadFluxDn(:)
                          ! 諠第弌螟ｧ豌励�蟆��縺ｫ逕ｱ譚･縺吶ｋ荳句髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2]
                          ! downward flux originated from emittance of planetary atmospheire [W m-2]
  real(8), allocatable :: r_RadNetFlux(:)
                          ! 諠第弌螟ｧ豌励�蟆��縺ｫ逕ｱ譚･縺吶ｋ豁｣蜻ｳ荳雁髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2]
                          ! net upward flux originated from emittance of planetary atmospheire [W m-2]
  real(8), allocatable :: r_SolFluxUp(:)
                          ! 荳ｭ蠢�弌蜈･蟆�↓逕ｱ譚･縺吶ｋ荳雁髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2]
                          ! upward flux originated from emittance of central star [W m-2]
  real(8), allocatable :: r_SolFluxDn(:)
                          ! 荳ｭ蠢�弌蜈･蟆�↓逕ｱ譚･縺吶ｋ荳句髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2]
                          ! downward flux originated from emittance of central star [W m-2]
  real(8), allocatable :: r_SolNetFlux(:)
                          ! 荳ｭ蠢�弌蜈･蟆�↓逕ｱ譚･縺吶ｋ豁｣蜻ｳ荳雁髄縺阪ヵ繝ｩ繝�け繧ｹ [W m-2]
                          ! net upward flux originated from emittance of central star [W m-2]

  ! namelist file 縺ｮ險ｭ螳�
  ! setting of namelist file
  !
  namelist /vflux_nml/ &
    & InFileName,      &
    & OutFileName

  ! namelist file 縺ｮ隱ｭ縺ｿ霎ｼ縺ｿ
  ! read namelist file
  !
  open(fnml, file='dcrtm.nml')
  read(fnml, nml=vflux_nml)
  close(fnml)
  !print *, MinWN, MaxWN, ResWN

  print *, InFileName
  print *, OutFileName

  ! 驟榊�繧ｵ繧､繧ｺ繧定ｪｿ縺ｹ繧�
  ! check array size
  !
  call HistoryGetPointer(''//InFileName//'', 'r_Press', a_Axis )

  kmax = size( a_Axis ) - 1

  Nullify( a_Axis )

  !print *, kmax

  call HistoryGetPointer(''//InFileName//'', 'w_WaveNum', a_Axis )

  wmax = size( a_Axis ) - 1

  ResWN = a_Axis(2) - a_Axis(1)

  Nullify( a_Axis )

  !print *, wmax
  !print *, ResWN

  !wmax = abs((MaxWN - MinWN)/ResWN)

  allocate( r_Press       (0:kmax        ) )
  allocate( rw_RadFluxUp  (0:kmax, 0:wmax) )
  allocate( rw_RadFluxDn  (0:kmax, 0:wmax) )
  allocate( rw_SolFluxUp  (0:kmax, 0:wmax) )
  allocate( rw_SolFluxDn  (0:kmax, 0:wmax) )
  allocate( r_RadFluxUp   (0:kmax        ) )
  allocate( r_RadFluxDn   (0:kmax        ) )
  allocate( r_SolFluxUp   (0:kmax        ) )
  allocate( r_SolFluxDn   (0:kmax        ) )
  allocate( r_RadNetFlux  (0:kmax        ) )
  allocate( r_SolNetFlux  (0:kmax        ) )

  ! 驩帷峩繝励Ο繝輔ぃ繧､繝ｫ縺ｮ隱ｭ縺ｿ霎ｼ縺ｿ
  ! read profile data
  !
  call HistoryGet(''//InFileName//'', 'r_Press', r_Press)
  call HistoryGet(''//InFileName//'', 'rw_RadFluxUp', rw_RadFluxUp)
  call HistoryGet(''//InFileName//'', 'rw_RadFluxDn', rw_RadFluxDn)
  call HistoryGet(''//InFileName//'', 'rw_SolFluxUp', rw_SolFluxUp)
  call HistoryGet(''//InFileName//'', 'rw_SolFluxDn', rw_SolFluxDn)

  ! 豕｢謨ｰ遨榊�
  ! summation flux
  !
  r_RadFluxUp(:) = (rw_RadFluxUp(:,0) + rw_RadFluxUp(:,wmax))*0.5_8
  r_RadFluxDn(:) = (rw_RadFluxDn(:,0) + rw_RadFluxDn(:,wmax))*0.5_8
  r_SolFluxUp(:) = (rw_SolFluxUp(:,0) + rw_SolFluxUp(:,wmax))*0.5_8
  r_SolFluxDn(:) = (rw_SolFluxDn(:,0) + rw_SolFluxDn(:,wmax))*0.5_8

  do k_w = 1, wmax-1
    r_RadFluxUp(:) = r_RadFluxUp(:) + rw_RadFluxUp(:,k_w)
    r_RadFluxDn(:) = r_RadFluxDn(:) + rw_RadFluxDn(:,k_w)
    r_SolFluxUp(:) = r_SolFluxUp(:) + rw_SolFluxUp(:,k_w)
    r_SolFluxDn(:) = r_SolFluxDn(:) + rw_SolFluxDn(:,k_w)
  end do

  r_RadFluxUp(:) = r_RadFluxUp(:)*ResWN
  r_RadFluxDn(:) = r_RadFluxDn(:)*ResWN
  r_SolFluxUp(:) = r_SolFluxUp(:)*ResWN
  r_SolFluxDn(:) = r_SolFluxDn(:)*ResWN

  r_RadNetFlux = r_RadFluxUp - r_RadFluxDn
  r_SolNetFlux = r_SolFluxUp - r_SolFluxDn

  ! 蜃ｺ蜉帙ヵ繧｡繧､繝ｫ縺ｮ險ｭ螳�
  ! setting of output file
  !
  call HistoryCreate(                                    &
    & file=''//OutFileName//'', &
    & title='vertical flux; '//OutFileName//'', &
    & source='sum_verticalprof.f90;'//InFileName//'', &
    & institution='DCRTM',       &
    & dims=(/'r_Press'/), dimsizes=(/kmax+1/),      &
    & longnames=(/'pressure'/),       &
    & units=(/'Pa'/)                    )

  call HistoryPut('r_Press',r_Press)

  call HistoryAddVariable( &
    & varname='r_RadFluxUp', dims=(/'r_Press'/), &
    & longname='upward flux(planet)', units='W m-2', xtype='double')
  call HistoryPut('r_RadFluxUp',r_RadFluxUp)

  call HistoryAddVariable( &
    & varname='r_RadFluxDn', dims=(/'r_Press'/), &
    & longname='upward flux(planet)', units='W m-2', xtype='double')
  call HistoryPut('r_RadFluxDn',r_RadFluxDn)

  call HistoryAddVariable( &
    & varname='r_SolFluxUp', dims=(/'r_Press'/), &
    & longname='upward flux(solar)', units='W m-2', xtype='double')
  call HistoryPut('r_SolFluxUp',r_SolFluxUp)

  call HistoryAddVariable( &
    & varname='r_SolFluxDn', dims=(/'r_Press'/), &
    & longname='upward flux(solar)', units='W m-2', xtype='double')
  call HistoryPut('r_SolFluxDn',r_SolFluxDn)

  call HistoryAddVariable( &
    & varname='r_RadNetFlux', dims=(/'r_Press'/), &
    & longname='net upward flux(planet)', units='W m-2', xtype='double')
  call HistoryPut('r_RadNetFlux',r_RadNetFlux)

  call HistoryAddVariable( &
    & varname='r_SolNetFlux', dims=(/'r_Press'/), &
    & longname='net upward flux(solar)', units='W m-2', xtype='double')
  call HistoryPut('r_SolNetFlux',r_SolNetFlux)

  call HistoryClose

end program sum_verticalflux