pro centroid_plot,n,brightm,faintm,brightc,faintc,brightp,faintp,help=help

if keyword_set(help) then begin
	print,'PURPOSE:'
	print,'	Makes arrays representing differences in centroids
	print,'	for brightest and faintest components'
	print,'HELP:'
	print,'	pro centroid_plot,n,brightm,faintm,brightc,faintc,brightp,faintp'
	print,'	inputs: n (number of sources,integer)
	print,'	outputs:brightm,faintm,brightc,faintc,brightp,faintp (double arrays)'
	print,'		(All are coordinates of component centroids)'
	print,'		Also produces Postscript files: centroid_mc.ps, centroid_mp.ps'
	print,'		and centroid_pc.ps'
	print,'AUTHOR:'
	print,'	Rick Pernak, Goddard Space Flight Center'
	return
endif

temp_mem = 'mem.csv';'/home/pernak/RHESSI/paper/spreadsheet/mem.csv'
temp_clean = 'clean.csv';'/home/pernak/RHESSI/paper/spreadsheet/clean.csv'
temp_pixon = 'pixon.csv';'/home/pernak/RHESSI/paper/spreadsheet/pixon.csv'

;restore,'/home/pernak/RHESSI/paper/figures/Clean/normal/idl_plots/clean_figs_legend.sav
;restore,'/home/pernak/RHESSI/paper/figures/MEM/mem_figs_legend.sav'
;restore,'/home/pernak/RHESSI/paper/figures/Pixon/pixon_figs_legend.sav'

openr,1,temp_mem
mem=strarr(n)
readf,1,mem
close,1

openr,1,temp_clean
clean=strarr(n)
readf,1,clean
close,1

openr,1,temp_pixon
pixon = strarr(n)
readf,1,pixon
close,1

;make MEM centroid arrays
for ctr=0,(n-1) do begin
	split = str_sep(mem(ctr),9B)
	xsource1 = double(split[5])-1 ;minus 1 arcsec b/c accomadating for .5 pixel shift in
	ysource1 = double(split[6])-1	;MEM maps
	xsource2 = double(split[7])-1
	ysource2 = double(split[8])-1
	if (ctr eq 0) then begin
		brightm_x = [xsource1]
		brightm_y = [ysource1]
		faintm_x = [xsource2]
		faintm_y = [ysource2]
	endif else begin
		brightm_x = [brightm_x,xsource1]
		brightm_y = [brightm_y,ysource1]
		faintm_x = [faintm_x,xsource2]
		faintm_y = [faintm_y,ysource2]
	endelse
endfor

;make bright and faint coordinate arrays for MEM
brightm = dblarr(2,n)
brightm[0,*] = brightm_x
brightm[1,*] = brightm_y
faintm = dblarr(2,n)
faintm[0,*] = faintm_x
faintm[1,*] = faintm_y

;make Clean centroid arrays
for ctr=0,(n-1) do begin
	split = str_sep(clean(ctr),9B)
	xsource1 = double(split[5])
	ysource1 = double(split[6])
	xsource2 = double(split[7])
	ysource2 = double(split[8])
	if (ctr eq 0) then begin
		brightc_x = [xsource1]
		brightc_y = [ysource1]
		faintc_x = [xsource2]
		faintc_y = [ysource2]
	endif else begin
		brightc_x = [brightc_x,xsource1]
		brightc_y = [brightc_y,ysource1]
		faintc_x = [faintc_x,xsource2]
		faintc_y = [faintc_y,ysource2]
	endelse
endfor

;make bright and faint coordinate arrays for Clean
brightc = dblarr(2,n)
brightc[0,*] = brightc_x
brightc[1,*] = brightc_y
faintc = dblarr(2,n)
faintc[0,*] = faintc_x
faintc[1,*] = faintc_y

;make Pixon centroid arrays
for ctr=0,(n-1) do begin
	split = str_sep(pixon(ctr),9B)
	;units: arcsec
	xsource1 = double(split[5])
	ysource1 = double(split[6])
	xsource2 = double(split[7])
	ysource2 = double(split[8])
	if (ctr eq 0) then begin
		brightp_x = [xsource1]
		brightp_y = [ysource1]
		faintp_x = [xsource2]
		faintp_y = [ysource2]
	endif else begin
		brightp_x = [brightp_x,xsource1]
		brightp_y = [brightp_y,ysource1]
		faintp_x = [faintp_x,xsource2]
		faintp_y = [faintp_y,ysource2]
	endelse
endfor

;make bright and faint coordinate arrays for Pixon
brightp = dblarr(2,n)
brightp[0,*] = brightp_x
brightp[1,*] = brightp_y
faintp = dblarr(2,n)
faintp[0,*] = faintp_x
faintp[1,*] = faintp_y

;MEM vs. Clean
dx1 = brightm_x - brightc_x
dy1 = brightm_y - brightc_y
dx2 = faintm_x - faintc_x
dy2 = faintm_y - faintc_y

xr=[-5,5]
yr=[-5,5] ;arbitrary x and y ranges

;make centroid difference plot
set_plot,'ps'
device,filename='centroid_mc.ps'
plot,dx1,dy1,$
	xr=xr, yr=yr, xtitle='Centroid dx (arcsec)',ytitle='Centroid dy (arcsec)',$
	psym=5,/iso,/xst,/yst,charsize=1.5
for ctr=0,n-1 do begin ;label with flare numbers
	xyouts,dx1[ctr],dy1[ctr],strcompress(ctr+1,/rem),charsize=1.5
endfor

oplot,dx2,dy2,psym=1
for ctr=0,n-1 do begin
	xyouts,dx2[ctr],dy2[ctr],strcompress(ctr+1,/rem),charsize=1.5
endfor

legend,['Brighter Component','Fainter Component'],psym=[5,1],box=0,charsize=1.5

;plot dashed lines through origin
oplot,[0,0],yr,linestyle=2
oplot,xr,[0,0],linestyle=2

;MEM vs. Clean RMS
r1 = dblarr(n)
r2 = dblarr(n)
for ctr=0,n-1 do begin
	r1[ctr] = dx1[ctr]^2 + dy1[ctr]^2 ;actually r^2 (distance squared),bright
	r2[ctr] = dx2[ctr]^2 + dy2[ctr]^2 ;weak components
endfor

;bright components' rms,creates string for plot legend
rms1 = sqrt(total(r1)/n)
rms1 = string(rms1,format='("RMS!Dbright!N =",1x,d8.6)')

;weak components' rms,ditto
rms2 = sqrt(total(r2)/n)
rms2 = string(rms2,format='("RMS!Dfaint!N =",1x,d8.6)')

legend,/bottom,[rms1,rms2],box=0,charsize=1.5
device,/close

;MEM vs. Pixon
dx1 = brightm_x - brightp_x
dy1 = brightm_y - brightp_y
dx2 = faintm_x - faintp_x
dy2 = faintm_y - faintp_y

;make centroid difference plot
set_plot,'ps'
device,filename='centroid_mp.ps'
plot,dx1,dy1,$
	xr=xr, yr=yr, xtitle='Centroid dx (arcsec)',ytitle='Centroid dy (arcsec)',$
	psym=5,/iso,/xst,/yst,charsize=1.5
for ctr=0,n-1 do begin ;label with flare numbers
	xyouts,dx1[ctr],dy1[ctr],strcompress(ctr+1,/rem),charsize=1.5
endfor

oplot,dx2,dy2,psym=1
for ctr=0,n-1 do begin
	xyouts,dx2[ctr],dy2[ctr],strcompress(ctr+1,/rem),charsize=1.5
endfor

legend,['Brighter Component','Fainter Component'],psym=[5,1],box=0,charsize=1.5
oplot,[0,0],yr,linestyle=2
oplot,xr,[0,0],linestyle=2

;MEM vs. Pixon RMS
r1 = dblarr(n)
r2 = dblarr(n)
for ctr=0,n-1 do begin
	r1[ctr] = dx1[ctr]^2 + dy1[ctr]^2 ;actually r^2 (distance squared),bright
	r2[ctr] = dx2[ctr]^2 + dy2[ctr]^2 ;weak components
endfor

;bright components' rms, string for plot legend
rms1 = sqrt(total(r1)/n)
rms1 = string(rms1,format='("RMS!Dbright!N =",1x,d9.6)')

;weak components' rms, string for plot legend
rms2 = sqrt(total(r2)/n)
rms2 = string(rms2,format='("RMS!Dfaint!N =",1x,d9.6)')

legend,/bottom,[rms1,rms2],box=0,charsize=1.5
device,/close

;Pixon vs. Clean
dx1 = brightp_x - brightc_x
dy1 = brightp_y - brightc_y
dx2 = faintp_x - faintc_x
dy2 = faintp_y - faintc_y

;make centroid difference plot
set_plot,'ps'
device,filename='centroid_pc.ps'
plot,dx1,dy1,$
	xr=xr, yr=yr, xtitle='Centroid dx (arcsec)',ytitle='Centroid dy (arcsec)',$
	psym=5,/iso,/xst,/yst,charsize=1.5
for ctr=0,n-1 do begin ;label with flare numbers
	xyouts,dx1[ctr],dy1[ctr],strcompress(ctr+1,/rem),charsize=1.5
endfor

oplot,dx2,dy2,psym=1
for ctr=0,n-1 do begin
	xyouts,dx2[ctr],dy2[ctr],strcompress(ctr+1,/rem),charsize=1.5
endfor

legend,['Brighter Component','Fainter Component'],psym=[5,1],box=0,charsize=1.5
oplot,[0,0],yr,linestyle=2
oplot,xr,[0,0],linestyle=2

;Pixon vs. Clean RMS
r1 = dblarr(n)
r2 = dblarr(n)
for ctr=0,n-1 do begin
	r1[ctr] = dx1[ctr]^2 + dy1[ctr]^2 ;actually r^2 (distance squared),bright
	r2[ctr] = dx2[ctr]^2 + dy2[ctr]^2 ;weak components
endfor

;bright components' rms,string for plot legend
rms1 = sqrt(total(r1)/n)
rms1 = string(rms1,format='("RMS!Dbright!N =",1x,d9.6)')

;weak components' rms,string for plot legend
rms2 = sqrt(total(r2)/n)
rms2 = string(rms2,format='("RMS!Dfaint!N =",1x,d9.6)')

legend,/bottom,[rms1,rms2],box=0,charsize=1.5
device,/close
set_plot,'x'

return
end