;
;HISTORY:
;	Written 29-May-96 by M.Morrison
;V1.1	30-May-96 (MDM) - Corrected the units - they were off by about 30%
;			- Added plotting of the orbit data
;V1.2	30-May-96 (MDM) - Added fitting the slope to the blueshift trend
;V1.3	12-Aug-96 (MDM) - Changed how the reading of the data was done
;			  (use rd_loi_mmad)
;V1.4	12-Aug-96 (MDM) - Removed slope/trend of velocity (makes FFT
;			  easier to interpret)
;V1.5	12-Aug-96 (MDM) - Removed marking ALT period
;V1.6	23-Jan-97 (MDM) - New RD_LOI_MMAD which has DOUBLE times (no slop)
;
progver = 'PLOT_LOI_96MIN.PRO  Ver 1.6
;
clearplot
if (n_elements(sttim) eq 0) then sttim = '1-jun'
if (n_elements(entim) eq 0) then entim = '10-jun'
if (n_elements(which) eq 0) then which = 'v'
wtit = 'LOI-' + strupcase(which)
if (which eq 'v') then wtit = wtit + ' (m/s)'
;
if (n_elements(timarr) eq 0) then rd_loi_mmad, which, sttim, entim, timarr, data0, clean=2
;
if (n_elements(data) eq 0) then begin
    coeff = poly_fit( int2secarr(timarr), data0, 1, yfit)
    data = data0 - yfit
end
;
!p.multi=[0,1,2]
utplot, timarr, data0, /ynozero, ytit='Average ' + wtit
outplot, timarr, yfit
utplot, timarr, data, /ynozero, ytit='Average '+wtit+' - Ramp Fit'
plottime & plottime, 0, 0, progver & pause
;
step = 60		;1 minute
;
plot_fft, timarr, data, ytit='Average '+wtit, step=step		;, mark_period=96*60.
pause
clearplot
;
;-------------------------------------------------------
;
fmt_timer, timarr, t1, t2
;
if (which eq 'v') then begin
  if (n_elements(scv) eq 0) then begin
    mat = rd_tfile('$MDI_CAL_INFO/soho_orbit.txt', 8, 2)
    timarr_scv = anytim2ints( fid2ex(mat(1,*)) )
    scv = float( mat(2,*) )
    ;
    x1 = int2secarr(timarr_scv, '1-jan-96')
    x2 = int2secarr(timarr, '1-jan-96')
    ;scv2 = spline(x1, scv, x2)
    scv3 = interpol(scv, x1, x2)
  end
  ;
  tit = 'LOI-V full average
  ytit2 = 'S/C Vel (m/sec)
  ytit3 = 'Avg Dop Vel - S/C Vel (m/sec)

  !p.multi=0
  utplot, timarr_scv, scv, tit='S/C Velocity', ytit=ytit2, psym=2
  outplot, timarr, scv3
  plottime & plottime, 0, 0, progver & pause

  coeff1 = poly_fit(  int2secarr(timarr)/86400., data0, 1, yfit1)
  coeff2 = poly_fit(  int2secarr(timarr)/86400., data0-scv3, 1, yfit2)

  !p.multi=[0,1,3]
  utplot, timarr, data0, /ynozero, ytit='Average LOI-V (m/s)', psym=3, timerange=[t1,t2]
  outplot, timarr, yfit1 & plottimes, .7, .8, 'Slope = ' + strtrim(coeff1(1),2) + ' m/sec per day'
  utplot, timarr, scv3, /ynozero, tit='S/C Velocity', ytit=ytit2, timerange=[t1,t2]
  utplot, timarr, data0-scv3, /ynozero, tit=tit, ytit=ytit3, psym=3, timerange=[t1,t2]
  outplot, timarr, yfit2 & plottimes, .7, .8, 'Slope = ' + strtrim(coeff2(1),2) + ' m/sec per day'
  plottime & plottime, 0, 0, progver & pause
  clearplot
end

;
;-------------------------------------------------------
;
if (n_elements(yy2) eq 0) then begin
    if (data_type(timarr) eq 8) then xx1 = int2secarr(timarr) $
                                else xx1 = timarr
    n1 = n_elements(timarr)
    yy1 = data
    dxx1 = deriv_arr(xx1)
    if (n_elements(step) eq 0) then step = min(dxx1)        ;sample step size
    print, 'Sample step size used is: ', step
    ;
    n2 = fix(max(xx1)/step)
    xx2 = findgen(n2)*step
    yy2 = interpol(yy1, xx1, xx2)

    ncycle = fix(n2/96.)

    yy2a = reform(yy2(0:96*ncycle-1),96,ncycle)
    ;;plot,total(yy2a,2)/ncycle,psym=10, /ynoz
end

if (n_elements(qokb) eq 0) then begin
    qok = intarr(n2)
    ss = tim2dset( anytim2ints(timarr(0), off=xx2), anytim2ints(timarr(0), off=xx1) )
    qok(ss) = 1
    qokb = reform(qok(0:96*ncycle-1),96,ncycle)
    totarr = fltarr(96)
    ntotarr = 0
    ;;for i=0,ncycle-1 do if (max(abs(yy2a(*,i)+2125)) gt 125) then qokb(0,i) = -99
    for i=0,ncycle-1 do if (total(qokb(*,i)) gt 70) then begin & ntotarr=ntotarr+1 & $
			totarr = totarr + yy2a(*,i) & end
    totarr = totarr/ntotarr
end
;
fmt_timer, timarr, t1, t2
plot, totarr, psym=10,/ynozero, tit=t1 + ' to ' + t2 + "  ' + wtit, $
			ytit='LOI Average of ' + strtrim(ntotarr,2) + ' 96 minute periods (m/sec)'

print, 'Do   qmarkalt=1    to see the 5 minute window marked (sometimes incorrectly)
if (keyword_set(qmarkalt)) then begin
    tarr = anytim2ex(t1)
    nmin = tarr(0)*60 + tarr(1)	;number of minutes into the day
    t_alt = nmin mod 96
    oplot, [1,1]*t_alt, !y.crange
    oplot, [1,1]*t_alt+6, !y.crange
end
plottime & plottime, 0, 0, progver & pause

print, 'Do    qdoall=1    to see each 96 minute individual plots'
!p.multi=[0,4,4]
iplot = 0
if (keyword_set(qdoall)) then for i=0,ncycle-1 do if (total(qokb(*,i)) gt 70) then begin
    plot, yy2a(*,i), /ynozero
    iplot = iplot+1
    if ((iplot mod 16) eq 0) then begin & plottime & plottime, 0, 0, progver & pause & end
end

print, 'Redo all:  delvar, timarr, data, yy2, qokb, scv

pprint
end