;+
;
; Project   : Densitycube, reconstruction, 3D, and s3drs software
;                   
; Name      : imageintersection
;               
; Purpose   : used to manipulate 3D densitycubes and do 3D reconstructions
;               
; Explanation: given two images, returns the sum of only their intersection regions
; you can normalize before adding by setting the 'normalize' flag
; (does not modify the original data)
;               
; Use       : standalone or via s3drs GUI
;    
; Inputs    : a, b (each a 2D data)
;               
; Outputs   : 
;
; Keywords  : normalize
;               
; Common    : none
;               
; Restrictions: none
;               
; Side effects: 
;               
; Category    : visualization
;               
; Prev. Hist. : None.
;
; Written     : Sandy Antunes, NRC, March-April 2009
;               
;-            
; given two images, returns the sum of only their intersection regions
; you can normalize before adding by setting the 'normalize' flag
; (does not modify the original data)

FUNCTION imageintersection,a,b,normalize=normalize

  normalize=KEYWORD_SET(normalize)

  N=(size(a))[1]
  aa=bytarr(N,N)
  bb=bytarr(N,N)
  cc=bytarr(N,N)
  aa=(a gt 0)
  bb=(b gt 0)
  cc=( (aa+bb) gt 1)
  c=fltarr(N,N)

  if (normalize) then begin
      am=max(a)
      bm=max(b)
      c=( (a/am) + (b/bm) )*cc
  end else begin
      c=(     a      +     b      )*cc
  end

  return,c

END

PRO test_intersection,normalize=normalize

  normalize=KEYWORD_SET(normalize)

  a=make_array(64,64,/float)
  b=make_array(64,64,/float)
  a(0:32,0:63)=0.
  a(32:48,0:63)=1000
  a(48:63,0:63)=2000

  b(0:63,0:32)=0.
  b(0:63,32:48)=5.
  b(0:63,48:63)=10.
;  a=a/max(a)
;  b=b/max(b)

  c=imageintersection(a,b,normalize=normalize)
  datum=[[[a]],[[b]],[[c]]]
  tv_multi,datum,title='intersection test'

END