Visibility workshop

;;CODE FOR DEMO:
 
	;RESTORE, 'visibility_bag_20021126.sav',/V
	RESTORE, 'visibility_bag_20021126_5min.sav',/V
	RESTORE, 'srctrack_20021126.sav',/V
 
vis00=vis0
;ss=WHERE(vis00.isc GE 7) & vis0=vis00[ss]
ss=WHERE(vis00.trange[0] GE anytim('2002/11/26 20:00') AND vis00.trange[0] LE anytim('2002/11/26 22:00')) & vis0=vis00[ss]
 
	nvis=N_ELEMENTS(vis0)
	phase_map_center=vis0[0].xyoffset	;; reference phase_map_center to fool software with
 
 
		.RUN
		vis1='bla'
		FOR i=0L, nvis-1 DO BEGIN
			IF i MOD 1000 EQ 0 THEN PRINT, i
			t=get_edges(vis0[i].trange,/MEAN)
			newx=INTERPOL(srctrack.X, srctrack.t, t)
			newy=INTERPOL(srctrack.Y, srctrack.t, t)
			tmp=hsi_vis_shift_mapcenter(vis0[i], [newx,newy])	;;changing the phase so that each visibility is centered on source's exact centroid position
			IF datatype(vis1) EQ 'STR' THEN vis1=tmp ELSE vis1=[vis1,tmp]
		ENDFOR;i
		END;.RUN
		vis1.xyoffset=phase_map_center	;; makes the s/w believe all the visibilities are referenced to the same xyoffset (phase_map_center)
 
	;RESTORE,'vis1.sav',/V
 
	!P.MULTI=[0,1,3]
	LOADCT,5
 
	ERRWEIGHT=0
 
	ss=WHERE(vis0.isc NE 6 AND vis0.isc GE 3) & vis0f=hsi_vis_combine(hsi_vis_edit(vis0[ss]), ERRWEIGHT=ERRWEIGHT)
	uv_smooth, vis0f, map, reconstructed_map_visibilities=rmv
	plot_map, /LIMB, map, /ISO
 
	ss=WHERE(vis1.isc NE 6 AND vis1.isc GE 3) & vis1f=hsi_vis_combine(hsi_vis_edit(vis1[ss]), ERRWEIGHT=ERRWEIGHT)
	uv_smooth, vis1f, map
	plot_map, /LIMB, map, /ISO
 
 
 
	;; compare vis0 (rchi2: 9.5) and vis1 (rchi2:5.8)
	hsi_vis_fwdfit, vis0f, CIRCLE=0, srcout=srcout, fitstddev=fitstddev

;Visibilities are a Data Compression Scheme
;
;1. Calib_eventlist is our first data compression - away from single events to time bins
;
;2. Stacked eventlist is a further compression - multiple rotations are epoch folded
;
;3. Visibilities are a further compression by a factor of 6 to 12, (half rotations)
;
;4. phase & amplitude fitting is inherently linear
;
;5. Inherent problem is tension between number of roll bins and minimum phase requirement
;	 in visibility determination.  phase gap -> throw away varying fractions of otherwise good data
;
;6. Once it's in visibilities you have hidden your systematic errors
;
;7. Viz normalization is the default.  We really need to study this
;
;8. For well separated sources with structure, the roll binning problem cannot be overcome for visibilities
;
;Let's play with visibilities
 
hessi_data, im=im
ov = im->get(/obj, class='hsi_visibility')
 
HELP, ov
;Get them all
 
v = ov->getdata(/dump)
HELP, v,/st
;Manipulate them
vs = v[0:99]
 
;Set them back
vnew = PTRARR(1)
vnew[0]= PTR_NEW(vs)
ov->setdata, vnew
 
;Look at what we have now
 
HELP, ov->getdata(/dump)

; Time range for your image of interest...
trange = '23-jul-02 00:31:'+['20','40']
 
; Energy ranges - number of intervals needs to be equal to number of
;   thermal components, and the fractional count contribution of each
;   component should be very different in each range...
erange = [[6.3,7.3], [17.,18.]]
 
; Fractional count contribution of ONE of the components (e.g. super-hot)
;   for each of the specified energy ranges above.  This needs to be
;   determined from your spectral fit (e.g. from OSPEX).  Fractional
;   contribution is based on COUNTS, not photons.
htpct = [0.604375, 0.948072]
 
; Determinant of the fractional contribution matrix (a square matrix
;   with dimensions equal to the number of thermal components, but which
;   is uniquely determined by only (N-1) x (N) numbers [since all rows
;   must sum to 1].
dn = htpct[1] - htpct[0]
 
; Set up the visibility object - you MUST turn off conjugation!
o = hsi_visibility(im_time_int=trange, im_energy_binning=erange[*,0], $
   xyoff=[-885,-230], det_index_mask=[0,0,1,1,1,1,0,1,1], vis_conj=0)
; Obtain visibilities for each energy range
v1 = o->getdata(energy_band=erange[*,0])
v2 = o->getdata(energy_band=erange[*,1])
 
; Normalize the visibilities to obtain "relative" visibilities
v1.obsvis /= v1.totflux  &  v2.obsvis /= v2.totflux
v1.sigamp /= v1.totflux  &  v2.sigamp /= v2.totflux
 
; Linearly combine raw (energy-space) visibilities using the inverse of
;   the fractional contribution matrix as weighting factors, to obtain
;   processed (temperature-space) visibilities for each thermal component.
vlo = hsi_vis_combine([hsi_vis_scale(v1, htpct[1]/dn), $
                       hsi_vis_scale(v2,-htpct[0]/dn)], /add,viscount=2)
vhi = hsi_vis_combine([hsi_vis_scale(v1,(htpct[1]-1)/dn), $
                       hsi_vis_scale(v2,(1-htpct[0])/dn)], /add,viscount=2)
 
; Plot the raw images and the processed images:
!p.MULTI=[0,2,1]
linecolors
uv_smooth, v1, m1  &  uv_smooth, v2, m2
uv_smooth, vlo, mlo  &  uv_smooth, vhi, mhi
plot_map, m1, /no_data, title='Raw Vis. [energy]'
  plot_map, m1, /over, /percent, levels=[50,75,95]
  plot_map, m2, /over, /percent, levels=[50,75,95], c_color=2
plot_map, mlo, /no_data, title='Thermal Vis. [temperature]'
  plot_map, mlo, /over, /percent, levels=[50,75,95]
  plot_map, mhi, /over, /percent, levels=[50,75,95], c_color=2
 
; END