!=======================================================================
!  エネルギー収支解析プログラム
!
!  2001/05/31  小高正嗣
!
!=======================================================================
      program energy
!-----------------------------------------------------------------------
      implicit real*8 (a-h, o-z)
!-----------------------------------------------------------------------
      parameter (ndim = 300)

      real*4 time(ndim), ETHERM(ndim), EKINET(ndim), ETOTAL(ndim),
     \     CTHKIN(ndim), DEKVIS(ndim), DEKNLV(ndim), EKPROG(ndim)

      real*4 hour, TMAX, ETMAX, ETMIN, EKMAX, EGMAX, EGMIN
      real*4 xvpmin, xvpmax, yvpmin, yvpmax

!-----------------------------------------------------------------------
      xvpmin = 0.2
      xvpmax = 0.8
      yvpmin = 0.15
      yvpmax = 0.6

      hour =3600.0

      open(11)
      read(11,*) idn
      close(11)

      open(10)

      do id = 1, idn
         read(10,*) time(id)
         read(10,*) ETHERM(id)
         read(10,*) EKINET(id)
         read(10,*) ETOTAL(id)
         read(10,*) CTHKIN(id)
         read(10,*) DEKVIS(id)
         read(10,*) DEKNLV(id)
      end do

      EKPROG(1) = CTHKIN(1) + DEKVIS(1) +DEKNLV(1)
      do id = 2, idn
         EKPROG(id) = EKPROG(id-1) + CTHKIN(id) + DEKVIS(id) +DEKNLV(id)
      end do

      close(10)

      do id = 1, idn
         time(id) = time(id) / hour
      end do

      ETMAX = 0. ; ETMIN = 0. ; EGMAX = 0. ; EGMIN = 0. ; EKMAX = 0. 

      do id = 1, idn
         ETMAX = max( ETMAX, ETOTAL(id) )
         ETMIN = min( ETMIN, ETOTAL(id) )
         EGMAX = max( EGMAX, CTHKIN(id), DEKVIS(id), DEKNLV(id) )
         EGMIN = min( EGMIN, CTHKIN(id), DEKVIS(id), DEKNLV(id) )
         EKMAX = max( EKMAX, EKINET(id) )
      end do

      call sgpwsn
      read(*,*) iws
      call gropn ( iws )
      call sldiv ( 'Y', 2, 2 )
      call sglset( 'LCORNER', .FALSE. )
      call sglset( 'LFULL', .TRUE. )

      call grfrm
      call grswnd( time(1), time(idn), ETMIN, ETMAX )
      call grsvpt( xvpmin, xvpmax, yvpmin, yvpmax )
      call grstrn( 1 )
      call grstrf
      call ussttl( 'Time', '[hour]', 'Energy', '[J]' ) 
      call usdaxs

      call sglset( 'LCHAR', .TRUE. )
      call sgspli( 24 )
      call sgsplc( 'total' )
      call sgplu ( idn, time, ETOTAL )
      call sgspli( 34 )
      call sgsplc( 'kinetic' )
      call sgplu ( idn, time, EKINET )
      call sgspli( 44 )
      call sgsplc( 'therm' )
      call sgplu ( idn, time, ETHERM )
      call sglset( 'LCHAR', .FALSE. )

      call grfrm
      call grswnd( time(1), time(idn), EGMIN, EGMAX )
      call grsvpt( xvpmin, xvpmax, yvpmin, yvpmax )
      call grstrn( 1 )
      call grstrf
      call ussttl( 'Time', '[hour]', 'Energy Generation', '[W]' ) 
      call usdaxs

      call sglset( 'LCHAR', .TRUE. )
      call sgspli( 24 )
      call sgsplc( 'T to K' )
      call sgplu ( idn, time, CTHKIN )
      call sgspli( 34 )
      call sgsplc( 'dis' )
      call sgplu ( idn, time, DEKVIS )
      call sgspli( 44 )
      call sgsplc( 'nlv' )
      call sgplu ( idn, time, DEKNLV )
      call sglset( 'LCHAR', .FALSE. )

      call grfrm
      call grswnd( time(1), time(idn), 0.0, EKMAX )
      call grsvpt( xvpmin, xvpmax, yvpmin, yvpmax )
      call grstrn( 1 )
      call grstrf
      call ussttl( 'Time', '[hour]', 'Kinetic Energy', '[W]' ) 
      call usdaxs

      call sglset( 'LCHAR', .TRUE. )
      call sgspli( 24 )
      call sgsplc( 'kinetic' )
      call sgplu ( idn, time, EKINET )
      call sgspli( 34 )
      call sgsplc( 'prog' )
      call sgplu ( idn, time, EKPROG )

      call sglset( 'LCHAR', .FALSE. )

      call grcls



      end