IGModel-SW 1.0

src/output_field_data.f90

説明を見る。
00001 
00011 module output_field_data
00012 
00013   ! モジュール引用 ; Use statements
00014   !
00015 
00016   ! 種類型パラメータ
00017   ! Kind type parameter 
00018   !
00019   use dc_types, only: DP, &    ! 倍精度実数型. Double precision. 
00020     &                 TOKEN, &
00021     &                 STRING !
00022 
00023   ! 座標変換
00024   ! Coordinate conversion
00025   !
00026   use igmcore_coordinate_conversion, only: &
00027     & orth_to_geo_pos, orth_to_geo_vec
00028 
00029   ! 正二十面体格子の管理
00030   ! Manager of icosahedral grid
00031   !
00032   use IcGrid2D_FVM_Manager, only: &
00033     & IcGrid2D_FVM, &
00034     & check_pole, get_EffSize_Min, get_EffSize_Max, &
00035     & RC_REGIONS_NUM, NOT_POLE_FLAG
00036 
00037   ! 正二十面体格子上の物理場の管理
00038   ! Manager of the physical field on icosahedral grid
00039   !
00040   use Field_IcGrid2D_Manager, only: &
00041     & Field_IcGrid2D, &
00042     & get_icgrid, &
00043     & Field_IcGrid2D_Init
00044   
00045   ! 正二十面体格子上の物理場データを NetCDF ファイルに書きだすクラス 
00046   ! Class to write the data of the physical fields on icosahedral grid into a NetCDF file
00047   !
00048   use IcGrid_ncWriter_mod, only: &
00049     & IcGrid_ncWriter, &
00050     & IcGrid_ncWriter_Init, open_ncFile, &
00051     & ncdef_simulation_parameter, ncdef_gridData, ncdef_FieldData, end_ncdef, &
00052     & write_GridData, write_FieldData, close_ncFile, &
00053     & increase_recorde_counter
00054 
00055   ! OpenMP
00056   !
00057   !
00058   use omp_lib
00059 
00060   ! 宣言文 ; Declaration statements
00061   !
00062   implicit none
00063   private
00064 
00065   ! 非公開変数
00066   ! Private variables
00067   !
00068 
00071   character(Token) :: ncFileName
00072   
00075   type(IcGrid_ncWriter), save :: writer
00076 
00079   type(IcGrid_ncWriter), save :: hs_writer
00080 
00083   type(Field_IcGrid2D), save :: UWind
00084   
00087   type(Field_IcGrid2D), save :: VWind
00088   
00091   type(Field_IcGrid2D), save :: WWind
00092 
00095   type(Field_IcGrid2D), save :: Height
00096 
00099   type(Field_IcGrid2D), save :: HeightS
00100 
00103   integer :: u_ncVarID
00104   
00107   integer :: v_ncVarID
00108   
00111   integer :: w_ncVarID
00112 
00115   integer :: h_ncVarID
00116 
00119   integer :: hs_ncVarID
00120 
00121   ! 公開関数
00122   ! Public procedures
00123   !
00124   public :: init_output_field_data, output_field, finalize_output_field_data
00125 
00126 contains
00127 
00128 !
00147 subroutine init_output_field_data( &
00148   & nc_filename, icgrid,           &            ! (in)
00149   & integration_time, time_step, output_tick )  ! (in)
00150 
00151   ! モジュール引用 ; Use statements
00152   !
00153   use dc_string, only: Replace
00154 
00155   ! 宣言文 ; Declaration statement
00156   !
00157   character(*), intent(in) :: nc_filename
00158   type(IcGrid2D_FVM), intent(in) :: icgrid
00159   real(DP), intent(in) :: integration_time
00160   real(DP), intent(in) :: time_step
00161   real(DP), intent(in) :: output_tick
00162 
00163   ! 作業変数
00164   ! Work variables
00165   !
00166   character(STRING) :: nc_hs_filename 
00167 
00168   ! 実行文 ; Executable statement
00169   !
00170   
00171   ncFileName = nc_filename
00172 
00173   ! IcGrid_ncWriter クラスのインスタンスを初期化. 
00174   call IcGrid_ncWriter_Init(writer, icgrid)
00175   call IcGrid_ncWriter_Init(hs_writer, icgrid)
00176 
00177   ! netCDF ファイルを開く. 
00178   call open_ncFile(writer, nc_filename)
00179   
00180   nc_hs_filename = trim(Replace(nc_filename, '.nc', '')) // '_hs.nc' 
00181   call open_ncFile(hs_writer, nc_hs_filename)
00182 
00183   ! netCDF ファイルへのアウトプットデータの定義
00184   !
00185 
00186   ! シミュレーションパラメータを NetCDF ファイルに書きだす. 
00187   call ncdef_Simulation_Parameter(writer, integration_time, time_step, output_tick)
00188   call ncdef_Simulation_Parameter(hs_writer, integration_time, time_step, output_tick)
00189 
00190   ! 格子点情報を NetCDF ファイルに定義する.  
00191   call ncdef_GridData(writer)
00192   call ncdef_GridData(hs_writer)
00193 
00194   ! 速度場と高度場の物理場データ情報を NetCDF ファイルに定義する.  
00195   call Field_IcGrid2D_Init(UWind, icgrid, 'u', 1, 'zonal_velocity', 'm/s')
00196   call Field_IcGrid2D_Init(VWind, icgrid, 'v', 1, 'meriodinal_velocity', 'm/s')
00197   call Field_IcGrid2D_Init(WWind, icgrid, 'w', 1, 'vertical_velocity', 'm/s')
00198   call Field_IcGrid2D_Init(Height, icgrid, 'h', 1, 'surface_height', 'm')
00199   call Field_IcGrid2D_Init(HeightS, icgrid, 'hs', 1, 'height_of_underlying_mountain', 'm')
00200 
00201   ! 各 Field_IcGrid2D クラスオブジェクトに対応付けられる NetCDF 内での変数管理 ID を取得する.
00202   u_ncVarID = ncdef_FieldData(writer, UWind)
00203   v_ncVarID = ncdef_FieldData(writer, VWind)
00204   w_ncVarID = ncdef_FieldData(writer, WWind)
00205   h_ncVarID = ncdef_FieldData(writer, Height)
00206   hs_ncVarID = ncdef_FieldData(hs_writer, HeightS)
00207 
00208   ! 物理場の定義を終了する. 
00209   call end_ncdef(writer)
00210   call end_ncdef(hs_writer)
00211 
00212   ! 格子情報を netCDF に書きこむ. 
00213   !
00214   call write_GridData(writer)
00215   call write_GridData(hs_writer)
00216 
00217 end subroutine init_output_field_data
00218 
00219 !
00236 subroutine output_field( &
00237   &  tstep, v_n, h_n, hs     &  ! (in)
00238   & )
00239 
00240   ! 宣言文 ; Declaration statement
00241   !
00242   integer, intent(in) :: tstep
00243   type(Field_IcGrid2D), intent(in) :: v_n
00244   type(Field_IcGrid2D), intent(in) :: h_n
00245   type(Field_IcGrid2D), intent(in), optional :: hs
00246 
00247   ! 作業変数 
00248   ! Work variables
00249   !
00250   type(IcGrid2D_FVM), pointer :: icgrid
00251   real(DP) :: s_v(3)
00252   integer :: rcID, i,j
00253   integer :: EMin, EMax
00254 
00255   ! 実行文 ; Executable statement
00256   !
00257 
00258   !write(*,*) 'ouput field: tstep=', tstep
00259   icgrid => get_icgrid(v_n)
00260 
00261   EMin = get_EffSize_Min(icgrid)
00262   EMax = get_EffSize_Max(icgrid)
00263 
00264   ! デカルト座標系で表されている速度場(ベクトル)データを, 地理座標系で表す. 
00265   do rcID=1, RC_REGIONS_NUM
00266     !$omp parallel do private(s_v)
00267     do j=EMin,EMax
00268       do i=EMin, EMax
00269         if ( check_pole(icgrid, rcID,i,j) == NOT_POLE_FLAG ) then
00270           s_v(:) = orth_to_geo_vec( v_n%val(rcID,i,j,:), orth_to_geo_pos(icgrid%rcs_AGrid(rcID,i,j,:)) )
00271           
00272           UWind%val(rcID,i,j,1) = s_v(1)
00273           VWind%val(rcID,i,j,1) = s_v(2)
00274           WWind%val(rcID,i,j,1) = s_v(3)
00275         end if
00276       end do
00277     end do
00278   end do
00279 
00280   Height%val(:,:,:,:) = h_n%val(:,:,:,:)
00281 
00282   ! 時間レベル n の速度場と高度場の情報を netCDF ファイルに書き出す. 
00283   call write_FieldData(writer, u_ncVarID, UWind)
00284   call write_FieldData(writer, v_ncVarID, VWind)
00285   call write_FieldData(writer, w_ncVarID, WWind)
00286   call write_FieldData(writer, h_ncVarID, Height)
00287 
00288   ! 地形の高度場を NetCDF ファイルに書きだす. 
00289   ! 地形の高度場については, 初期場に対してのみ書きだす. 
00290   if( tstep == 0 .and. present(hs) ) then
00291     HeightS%val(:,:,:,:) = hs%val(:,:,:,:)
00292     call write_FieldData(hs_writer, hs_ncVarID, HeightS)
00293     call close_ncFile(hs_writer)
00294   end if
00295 
00296   ! 時間レベル記録カウンタをインクリメントする. 
00297   call increase_recorde_counter(writer)
00298 
00299 end subroutine output_field
00300 
00301 !
00307 subroutine finalize_output_field_data()
00308 
00309   ! 実行文 ; Executable statement
00310   !
00311   
00312   ! NetCDF ファイルの書き出しクラスのオブジェクトを最終化する. 
00313   call close_ncFile(writer)
00314 
00315 end subroutine finalize_output_field_data
00316 
00317 end module output_field_data
 全て クラス ネームスペース ファイル 関数 変数