!
! 各矩形領域の物理量データを, vtk 形式でファイルに書き出すためのサブルーチンを提供するモジュール
!
module rcregion_vtkwriter
  implicit none
  private

  public :: init_rcregion_vtkwriter, close_rcregion_vtkwriter
  public :: begin_rcregion_grid_data, write_rcregion_grid_data
  public :: begin_point_data, write_point_data_scalars, write_point_data_vectors

  integer, save :: deviceID                       !

  ! vtk フォーマットが定める定数
  integer, parameter :: VTK_VERTEX = 1
  integer, parameter :: VTK_TRIANGLE_STRIP = 6

contains
subroutine init_rcregion_vtkwriter(device_id, rcID, dir_path, identifier )
  !
  ! rcregion_vtkwriter モジュールを初期化する.
  ! 書き出すファイルを作成し, vtk 形式のファイルに必要なヘッダーをまず書き出す.
  !
  integer, intent(in) :: device_id       ! 装置番号
  integer, intent(in) :: rcID            ! 対象としている矩形領域の id
  character(*), intent(in) :: dir_path   ! データファイルを書き出すディレクトリのパス
  character(*), intent(in) :: identifier ! 同矩形領域のデータファイル名を区別するための名前(例えばタイムレベル)

  character(len=2) c_rcID

  deviceID = device_id
  write(c_rcID, '(i2.2)') rcID

  open(deviceID, FILE=dir_path // 'field' // c_rcID // '-' // identifier // '.vtk' )
  write(deviceID,'(a)') '# vtk DataFile Version 2.0'
  write(deviceID,'(a)') 'IGModel Output Data'
  write(deviceID,'(a)') 'ASCII'
end subroutine init_rcregion_vtkwriter


subroutine begin_rcregion_grid_data(dataNum)
  !
  ! 格子データの書き出し始める前に呼び出すサブルーチン.
  ! 格子データの前に必要なヘッダーを書き出す.
  !
  integer, intent(in) :: dataNum !　対象としている矩形領域の格子データ数

  ! データセットが, 非構造格子であることを指定
  write(deviceID,'(a)') 'DATASET UNSTRUCTURED_GRID'
  ! データ数(格子点数)を指定
  write(deviceID,'(a,i5,a6)') 'POINTS', dataNum, ' float'

end subroutine begin_rcregion_grid_data


subroutine write_rcregion_grid_data(std_grid)
  !
  ! 格子データを書き出す.
  !
  real(8), intent(in) :: std_grid(:,:,:) ! 対象としている矩形領域の格子データを保持している配列

  integer i, j, imax, jmax
  integer cell_num, cell_pt_num

  imax = size(std_grid, 1)
  jmax = size(std_grid, 2)

  do i=1, imax
    do j=1, jmax
      write(deviceID,'(3e15.5e3)') std_grid(i, j, 1:3)
    end do
  end do
  write(deviceID, *)

  ! 要素の形状(セル)の設定を書き出す.
  cell_num = jmax-1
  cell_pt_num =2*imax
  write(deviceID, '(a, i5, i5)') 'CELLS', cell_num, (cell_pt_num+1)*cell_num
  do j=0, cell_num-1
    write(deviceID, '(i5)', advance='no') cell_pt_num
    do i=0, imax-1
      write(deviceID, '(2i5)', advance='no') i*jmax+j, i*jmax+j+1
    end do
    write(deviceID,*)
  end do
  write(deviceID, *)

  write(deviceID, '(a, i5)') 'CELL_TYPES', cell_num
  do i=0, cell_num-1
    write(deviceID, '(i1)') VTK_TRIANGLE_STRIP
  end do
  write(deviceID, *)

end subroutine write_rcregion_grid_data

subroutine begin_point_data(dataNum)
  !
  ! 格子点上のデータを書き出す前に呼び出すサブルーチン
  !
  integer, intent(in) :: dataNum  ! 書き込むデータ数

  write(deviceID, '(a, i5)') 'POINT_DATA', dataNum
end subroutine begin_point_data

subroutine write_point_data_scalars(dataName, scalar_data)
  !
  ! 格子点上のスカラーデータを書き出す
  !
  character(*), intent(in) :: dataName
  real(8), intent(in) :: scalar_data(:, :)

  integer i, j, imax, jmax

  imax = size(scalar_data, 1)
  jmax = size(scalar_data, 2)


  write(deviceID, '(a,a,a,i2)') 'SCALARS ', dataName, ' float ', VTK_VERTEX
  write(deviceID, '(a)') 'LOOKUP_TABLE defalut'

  do i=1, imax
    do j=1, jmax
      write(deviceID, '(e15.5e3)') scalar_data(i,j)
    end do
  end do
  write(deviceID, *)

end subroutine write_point_data_scalars

subroutine write_point_data_vectors( dataName, vector_data )
  !
  ! 格子点上のベクトルデータを書き出す
  !
  character(*), intent(in) :: dataName
  real(8), intent(in) :: vector_data(:, :, :)

  integer i, j, imax, jmax

  imax = size(vector_data, 1)
  jmax = size(vector_data, 2)
  write(deviceID, '(a,a,a)') 'VECTORS ', dataName, ' float '

  do i=1, imax
    do j=1, jmax
      write(deviceID, '(3e15.5e3)') vector_data(i,j, 1:3)
    end do
  end do
  write(deviceID, *)

end subroutine write_point_data_vectors

subroutine close_rcregion_vtkwriter()
  !
  ! ファイルを閉じる.
  !

  close(deviceID)
end subroutine close_rcregion_vtkwriter

end module rcregion_vtkwriter
