; Plot U from GOES data

; Initialize - 0 for windows, 1 for postscript (portrait), 2 for postscript (landscape)
; type ps=0, 1, or 2 on command line
checkvar,ps,0
initialize, ps

restore,'idlsave_goes031103.dat'

; Calculate start time

; Time range to plot
s = 40. & e = n_elements(tarray) - 1000
stime = anytim(utbase+tarray(s),/atime)
etime = anytim(utbase+tarray(e),/atime)

;-----------------------------------------------------------

c = rd_ascii('flux_output_nov3_12-25_d345678_c50_changed.txt')
c = c[1:*]
b = str2cols(c)
box = fix(reform(b[0,*]))
areas = float(reform(b[5,*]))
q = where(box eq 0)
b0 = b[*,q]
box = (box[q])
areas = (areas[q])
box = fix(reform(b0[0,*]))
areas = float(reform(b0[5,*]))
statimes = reform(b0[1,*] + ' ' + b0[2,*])
endtimes=reform(b0[1,*] + ' ' + b0[3,*])
stime_sec = anytim(reform(statimes))
etime_sec = anytim(reform(endtimes))

goestimes = tarray+utbase
vec = goestimes
vol = (4./3.) * (1. / sqrt (!pi))
im_vol = vol * areas^(1.5)*72498958.3^3

v = fltarr(n_elements(tarray)) + 4.7e27

for i=0, n_elements(areas)-1 do begin
	i1 = value_locate(vec, stime_sec[i])
	i2 = value_locate(vec, etime_sec[i])
;	print, i1,i2, im_vol[i]
	 v[i1:i2] = im_vol[i]

endfor


;----------------------------------------------------------


;Plot thermal energy U = 3 n k T V erg
; nV = SQRT(EM * V)
; Boltzmann constant - k = 1.38e-16 erg K^-1
   k = 1.38e-16
; i.e., U = 4.14e-16 T(K) sqrt(EM(cm^-3) * V(cm^-3)) ergs
; U = 4.14e-16 tempr * 1.e6 (K) sqrt(emis*1.e49(cm^-3) * V(cm^3)) ergs
; U = 4.14 e-10 * sqrt(10) * 1.e24 *tempr * sqrt(emis * V) ergs
; U = 1.309e15 * tempr * sqrt(emis * V) erg
; V = source volume in cm^3
;reads ascii file and returns it as a text array

   U = 1.309e15 * tempr * sqrt(emis * V)


utplot,tarray[s:e],U[s:e],utbase, $
   /ylog,ytitle='Thermal Energy U in ergs'
;   color=0,background=255,xticklen = 1,$
;   xtitle = ' ', ymargin = [0,0], yrange=[1.e28,1.e32], xmargin = [20,0]

if ps gt 0 then device,/close
set_plot,'win'

end