GOTO, Donereading

num = 238
gb = intarr(1237, 992, num)
bc = intarr(1237, 992, num)
rc = intarr(1237, 992, num)
sall = intarr(2, num)
nsum = 3
sh = intarr(2, num)
scl = 5000

fg = findfile('../FITS/GBand/AR*')
fb = findfile('../FITS/BlueCont/AR*')
fr = findfile('../FITS/RedCont/AR*')

window, 1, xs=1237, ys=992

FOR i=0, num-1 DO BEGIN
      ;;; read G band
    rdfits, p, fg(i)
    p1 = reform(p(*, *, 0))
    sgb = size(p1, /dim)
    gb(0, 0, i) = fix(scl*p1)
      ;;; read blue cont
    rdfits, p, fb(i)
    p1 = rot(reform(p(*, *, 0)), 0.95399, 1.00220, cub=-.5)
    sbc = size(p1, /dim)
    bc(0, 0, i) = fix(scl*p1)
      ;;; read red cont
    rdfits, p, fr(i)
    p1 = rot(reform(p(*, *, 0)), 1.13294, 1.00427, cub = -.5)
    src = size(p1, /dim)
    rc(0, 0, i) = fix(scl*p1)
    
      ;;; now shift all triples to match and cut the size to be the same
    
    ssh = intarr(2, 3)
    ssh(*, 1) = shc(gb(*, *, i), bc(*, *, i), /n2, /filt)
    ssh(*, 2) = shc(gb(*, *, i), rc(*, *, i), /n2, /filt)
    ssh(0, *) = ssh(0, *)-max(ssh(0, *))
    ssh(1, *) = ssh(1, *)-max(ssh(1, *))
    gb(*, *, i) = shift(gb(*, *, i), ssh(0, 0), ssh(1, 0))
    bc(*, *, i) = shift(bc(*, *, i), ssh(0, 1), ssh(1, 1))
    rc(*, *, i) = shift(rc(*, *, i), ssh(0, 2), ssh(1, 2))
    sall(*, i) = (sgb+ssh(*, 0)) < (sbc+ssh(*, 1)) < (src+ssh(*, 2))
    
    tvscl, gb(*, *, i) & wait, .5
    tvscl, bc(*, *, i) & wait, .5
    tvscl, rc(*, *, i) & empty
    
    IF i EQ 0 THEN BEGIN
        pp = fltarr(nsum, sall(0, 0), sall(1, 0))
        pp(0, *, *) = float(gb(0:sall(0, 0)-1, 0:sall(1, 0)-1, 0))+ $
          bc(0:sall(0, 0)-1, 0:sall(1, 0)-1, 0)+ $
          rc(0:sall(0, 0)-1, 0:sall(1, 0)-1, 0)
    ENDIF ELSE BEGIN
        p1=float(gb(0:sall(0, i)-1, 0:sall(1, i)-1, i))+ $
          bc(0:sall(0, i)-1, 0:sall(1, i)-1, i)+ $
          rc(0:sall(0, i)-1, 0:sall(1, i)-1, i)
        sh(*, i) = shc(total(pp, 1), p1, /filt, /N2)
        sx = (sall(0, i) < sall(0, 0))
        sy = (sall(1, i) < sall(1, 0))
        p1 = p1(0:sx-1, 0:sy-1)
        pp(i MOD nsum, 0:sx-1, 0:sy-1) = shift(p1, sh(0, i), sh(1, i))
        print, i, sh(*, i)
    ENDELSE
    
ENDFOR

sx = max(sh(0, *))-sh(0, *)
sy = max(sh(1, *))-sh(1, *)

resx = min(sall(0, *)-sx)
resy = min(sall(1, *)-sy)

  ;;; mave even dimensions
resx = resx/2*2
resy = resy/2*2


FOR i=0, num-1 DO BEGIN
    gb(*, *, i) = shift(gb(*, *, i), -sx(i), -sy(i))
    bc(*, *, i) = shift(bc(*, *, i), -sx(i), -sy(i))
    rc(*, *, i) = shift(rc(*, *, i), -sx(i), -sy(i))
ENDFOR

gb = (temporary(gb))(0:resx-1, 0:resy-1, *)
bc = (temporary(bc))(0:resx-1, 0:resy-1, *)
rc = (temporary(rc))(0:resx-1, 0:resy-1, *)

enhance, gb
enhance, bc
enhance, rc

openw, 1, 'gb_enh_1118x868x238' & writeu, 1, gb & close, 1
openw, 1, 'bc_enh_1118x868x238' & writeu, 1, bc & close, 1
openw, 1, 'rc_enh_1118x868x238' & writeu, 1, rc & close, 1

openr, 1, 'gb_enh_1118x868x238' & readu, 1, gb & close, 1
openr, 1, 'bc_enh_1118x868x238' & readu, 1, bc & close, 1
openr, 1, 'rc_enh_1118x868x238' & readu, 1, rc & close, 1


  ;;; now the rubbersheeting
FOR i=1, num-2 DO BEGIN
    ref = total(gb(*, *, i-1:i+1), 3) + total(bc(*, *, i-1:i+1), 3) + $
      total(rc(*, *, i-1:i+1), 3)
    gb(*, *, i)=remap_image(ref, gb(*, *, i), 24, 10, $
                            plot=3, range=3, /same, /quint)
    bc(*, *, i)=remap_image(ref, bc(*, *, i), 24, 10, $
                            plot=3, range=3, /same, /quint)
    rc(*, *, i)=remap_image(ref, rc(*, *, i), 24, 10, $
                            plot=3, range=3, /same, /quint)

ENDFOR

openw, 1, 'gb_rub_1118x868x238' & writeu, 1, gb & close, 1
openw, 1, 'bc_rub_1118x868x238' & writeu, 1, bc & close, 1
openw, 1, 'rc_rub_1118x868x238' & writeu, 1, rc & close, 1

Donereading:
mkmov, ca, '../FITS/CaIIH/AR*'
openw, 1, 'ca_1148x876x238' & writeu, 1, ca & close, 1
enhance, ca
openw, 1, 'ca_enh_1148x876x238' & writeu, 1, ca & close, 1
rubb, ca, cell=24, step=10, plot=3, /same, range=3, /quint
FOR i=0, num-1 DO $
  ca(*, *, i) = 5000.*ca(*, *, i)/avg(ca(20:1127, 20:855, i))+.5
openw, 1, 'ca_rub_1148x876x238' & writeu, 1, ca & close, 1

pic = histo_opt(total(ca(*, *, 135:137), 3))
ref_pic = total(gb(*, *, 135:137), 3)
; .run find_points
;  618  482   0.0185966   1.0028214 -24.16 -21.09
;
; rotation is in radian!  This is 1.06 degrees

FOR i=0, num-1 DO ca(*, *, i)=$
  rot(ca(*, *, i), 1.06, 1./0.996585, 1148/2+24.16, 876/2+21.09, cub=.5)

ca = (temporary(ca))(0:1117, 0:853, *)
gb = (temporary(gb))(*, 0:853, *)
bc = (temporary(bc))(*, 0:853, *)
rc = (temporary(rc))(*, 0:853, *)

openw, 1, 'AR10407-gb-0000-20030714-cube-1118x854x238-dt30'
writeu, 1, gb & close, 1

openw, 1, 'AR10407-ca-0000-20030714-cube-1118x854x238-dt30'
writeu, 1, ca & close, 1

openw, 1, 'AR10407-bc-0000-20030714-cube-1118x854x238-dt30'
writeu, 1, bc & close, 1

openw, 1, 'AR10407-rc-0000-20030714-cube-1118x854x238-dt30'
writeu, 1, rc & close, 1


conefilt_file, 'AR10407-gb-0000-20030714-cube-1118x854x238-dt30', $
  1118, 854, 238, 51.475, 30, 7, partsize=[576, 512], $
  resfile='AR10407-gb-0000-20030714-cube-1118x854x238-dt30-lo7'
conefilt_file, 'AR10407-bc-0000-20030714-cube-1118x854x238-dt30', $
  1118, 854, 238, 51.475, 30, 7, partsize=[576, 512], $
  resfile='AR10407-bc-0000-20030714-cube-1118x854x238-dt30-lo7'
conefilt_file, 'AR10407-rc-0000-20030714-cube-1118x854x238-dt30', $
  1118, 854, 238, 51.475, 30, 7, partsize=[576, 512], $
  resfile='AR10407-rc-0000-20030714-cube-1118x854x238-dt30-lo7'
conefilt_file, 'AR10407-ca-0000-20030714-cube-1118x854x238-dt30', $
  1118, 854, 238, 51.475, 30, 7, partsize=[576, 512], $
  resfile='AR10407-ca-0000-20030714-cube-1118x854x238-dt30-lo7'


tim = lindgen(238)*30
tim0 = 3600l*8+60l*22+30
tt = tim+tim0
hr = tt/3600
mn = (tt MOD 3600)/60
sc = tt MOD 60
tt = nnumber(hr, 2)+":"+nnumber(mn, 2)+":"+nnumber(sc, 2)

FOR i=0, num-1 DO BEGIN
    tv, m1(*, *, i)
    clock, tt(i), /dev, siz=76, pos=[5, 5]
    m1(*, *, i) = tvrd()
ENDFOR

mm = intarr(1118, 854, 238)
m1 = bytarr(1024, 768, 238)
openr, 1, 'AR10407-rc-0000-20030714-cube-1118x854x238-dt30-lo7'
readu, 1, mm & close, 1
mm1 = scale_array(rebin(mm, 1118/2, 854/2, 238))
m1(512:*, 0:383, *) = mm1(30:541, 10:393, *)
openr, 1, 'AR10407-bc-0000-20030714-cube-1118x854x238-dt30-lo7'
readu, 1, mm & close, 1
mm1 = scale_array(rebin(mm, 1118/2, 854/2, 238))
m1(512:*, 384:*, *) = mm1(30:541, 10:393, *)
openr, 1, 'AR10407-gb-0000-20030714-cube-1118x854x238-dt30-lo7'
readu, 1, mm & close, 1
mm1 = scale_array(rebin(mm, 1118/2, 854/2, 238))
m1(0:511, 384:*, *) = mm1(30:541, 10:393, *)
openr, 1, 'AR10407-ca-0000-20030714-cube-1118x854x238-dt30-lo7'
readu, 1, mm & close, 1
mm1 = scale_array(rebin(mm, 1118/2, 854/2, 238))
m1(0:511, 0:383, *) = mm1(30:541, 10:393, *)


FOR i=0, num-1 DO BEGIN
    tv, m1(*, *, i)
    clock, tt(i), /dev, siz=76, pos=[512, 384]-38
    xyouts, 10, 10, /dev, 'Ca II H', charsiz=2, charthick=2
    xyouts, 1013, 10, /dev, 'Red continuum', charsiz=2, charthick=2, alig=1
    xyouts, 10, 740, /dev, 'G band', charsiz=2, charthick=2
    xyouts, 1013, 740, /dev, 'Blue continuum', charsiz=2, charthick=2, alig=1
    m1(*, *, i) = tvrd()
ENDFOR

write_mpeg,'AR10407-all-20030714-lo7.mpg',m1,/high,tmp='/data/pit/tmp/mpegtmp'


END