pro perp_dist,help=help

if keyword_set(help) then begin
	print,'PURPOSE:'
	print,'	find perpendicular distance from y=x line for flux vs. flux plots '
	print,'	then calculates rms for the three plots and sigma (distance/rms) '
	print,'	for 1 flare'
	print,'HELP:'
	print,'	pro perp_dist'
	print,'	no inputs, output is printed on screen'
	print,'AUTHOR:'
	print,'	Rick Pernak, Goddard Space Flight Center'
	return
endif


openr,1,'/home/pernak/RHESSI/paper/spreadsheet/mem.csv'
openr,2,'/home/pernak/RHESSI/paper/spreadsheet/clean.csv'
openr,3,'/home/pernak/RHESSI/paper/spreadsheet/pixon.csv'

mem_arr = strarr(25)
clean_arr = strarr(25)
pixon_arr = strarr(25)

readf,1,mem_arr
readf,2,clean_arr
readf,3,pixon_arr

close,1,2,3

;combine all arrays for looping, "method" array
method = strarr(3,25)
method[0,*] = mem_arr
method[1,*] = clean_arr
method[2,*] = pixon_arr

;initialize flux arrays, m:mem, c:clean, p:pixon, 1:bright, 2:faint
m1 = dblarr(25)
m2 = dblarr(25)
c1 = dblarr(25)
c2 = dblarr(25)
p1 = dblarr(25)
p2 = dblarr(25)

;distance arrays
mc1 = dblarr(25)
mc2 = dblarr(25)
mp1 = dblarr(25)
mp2 = dblarr(25)
pc1 = dblarr(25)
pc2 = dblarr(25)

for n=0,2 do begin ;loop over method array
	for ctr=0,24 do begin ;loop over sources
		temp = method[n,*]
		split=str_sep(temp(ctr),'	')
		bright = split[1]
		faint = split[2]
		if (n eq 0) then begin
			m1[ctr] = bright
			m2[ctr] = faint
		endif else if (n eq 1) then begin
			c1[ctr] = bright
			c2[ctr] = faint
		endif else begin
			p1[ctr] = bright
			p2[ctr] = faint
		endelse
	endfor
endfor

;distance of flare on MEM vs. Clean flux from y=x
mc1 = abs(alog10(m1)-alog10(c1))/sqrt(2)
mc2 = abs(alog10(m2)-alog10(c2))/sqrt(2)
rms_mc1 = sqrt(mean(mc1)^2)
rms_mc2 = sqrt(mean(mc2)^2)

;distance of flare on MEM vs. Pixon flux from y=x
mp1 = abs(alog10(m1)-alog10(p1))/sqrt(2)
mp2 = abs(alog10(m2)-alog10(p2))/sqrt(2)
rms_mp1 = sqrt(mean(mp1)^2)
rms_mp2 = sqrt(mean(mp2)^2)

;distance of flare on Pixon vs. Clean flux from y=x
pc1 = abs(alog10(p1)-alog10(c1))/sqrt(2)
pc2 = abs(alog10(p2)-alog10(c2))/sqrt(2)
rms_pc1 = sqrt(mean(pc1)^2)
rms_pc2 = sqrt(mean(pc2)^2)

;in "log10 flux space"
print,'RMS mc1: ',rms_mc1
print,'RMS mc2: ',rms_mc2
print,'RMS mc average: ',(rms_mc1+rms_mc2)/2
print,'RMS mp1: ',rms_mp1
print,'RMS mp2: ',rms_mp2
print,'RMS mp average: ',(rms_mp1+rms_mp2)/2
print,'RMS pc1: ',rms_pc1
print,'RMS pc2: ',rms_pc2
print,'RMS pc average: ',(rms_pc1+rms_pc2)/2

print,'Flare 19 mp1/rms: ',mp1[18]/rms_mp1
print,'Flare 19 mp2/rms: ',mp2[18]/rms_mp2
print,'Flare 19 pc1/rms: ',pc1[18]/rms_pc1
print,'Flare 19 pc2/rms: ',pc2[18]/rms_pc2

return
end