!CCXで計算したトルク−ねじれ角、トルク-滑り角 関係をgnuplot用にまとめる

! 実行例：
! $ ./tth<toruku.dat>tth.d
! gnuplot< pl 'tth.d' using 1:4 w lp,'tth.d' using 2:4 w lp

implicit none
!---入力---
integer :: nstep=100   !     :出力ステップ数(CGXで確認)
integer :: nz=1000    !     :PC鋼棒の軸方向の分割数
real    :: tz=230.d-6 ! MNm :導入トルク
!---
real(8),allocatable :: thsube(:) ! rad :PC鋼棒の固定側の先端の回転角
real(8),allocatable :: thsai(:)  ! rad :PC鋼棒の載荷点の回転角
real(8),allocatable :: thnezi(:) ! rad :PC鋼棒のねじり角
real(8),allocatable :: t(:)      ! MNm :荷重
integer    i,j,icut
real       rcut
character  ccut*20
real(8)    step,pi
allocate(thsube(nstep),thsai(nstep),thnezi(nstep),t(nstep))

pi=2.d0*asin(1.d0)
thsube=0.d0
thsai=0.d0
thnezi=0.d0
t=0.d0
print*,'# th_nezi(degree) th_suberi(*100 degree) th_sai(degree) T(N)'!出力ファイルのヘッダー
do j=1,nstep
   read(*,'()')!datファイルのデータの読み飛ばし
   read(*,'()')!datファイルのデータの読み飛ばし
   read(*,'()')!datファイルのデータの読み飛ばし
   read(*,*) ccut,icut,rcut,rcut,thsube(j)!PC鋼棒の滑り角の読み込み
   do i=1,9
      read(*,'()')!datファイルのデータの読み飛ばし
   enddo
   !↓のdatファイルは、初めの文字と整数がくっついていることがあるので、read文は要注意!!
   read(*,*) ccut,rcut,rcut,thsai(j)!PC鋼棒の先端のねじれ角の読み込み
   !print'(1pe13.5,1pe13.5)',thsai(j),thsube(j)
   read(*,'()')!datファイルのデータの読み飛ばし
   read(*,'()')!datファイルのデータの読み飛ばし
   !datファイルの荷重の割合を読み込み
   read(*,*) ccut,ccut,ccut,ccut,ccut,ccut,ccut,ccut,ccut,step
   t(j)=(step-1.d0)*tz
   thnezi(j)=thsai(j)-thsube(j)
   do i=1,6
      read(*,'()')
   enddo
   !gnuplot用のデータの出力
   print'(1pe13.5,1pe13.5,1pe13.5,1pe13.5)',(thnezi(j)-thnezi(1))*180.d0/pi,(thsube(j)-thsube(1))*180.d2/pi,&
       &(thsai(j)-thsai(1))*180.d0/pi,t(j)*1.d6
enddo

end
