; csg.pro - Defines a surface with Constructive Solid Geometry (CSG) from
;           other surfaces
;
; ----------------------------------------------------------------------
; CONSTRUCTOR
; ----------------------------------------------------------------------
;Inputs:
;  elements - array of object references to other surfaces.
;  is_union - boolean scalar. true if the CSG is the union of the child
;             objects, false if the CSG is the intersection of the
;             child objects
function csg::init,elements,is_union
  ; init superclass
  status = self->surface::init()
  self.elements=ptr_new(elements);
  self.is_union=is_union;
  return, status
end

;function is_inside
;Inputs:
;  R - a vector to test. May be a grid
;Returns:
;  1 if given vector is inside the surface, 0 if not. In general, this should
;  return true if the point is exactly on the surface.
function csg::is_inside,r
  ;If the exception element is negative, it will never match and
  ;this will perform an is_inside on all elements
  return, self->is_inside_except(r,-1)
end

;function is_inside_except
;Inputs:
;  R - a vector to test. May be a grid
;  except - index of the element array to NOT check for insideness.
;Returns:
;  1 if given vector is inside the surface, 0 if not, with the exception above.
;
;  Generally this is used with intersection checking. The calculated intersection
;  point is very near the surface of the element which is intersected, but may be
;  inside, on, or outside, depending on floating point round-off, so we should not
;  check that element. Besides, we already know that this point is (should be)
;  exactly on the surface.
function csg::is_inside_except,r,except
  result=bytarr(shape_grid(r))+(~self.is_union)
  for i=0,n_elements((*(self.elements)))-1 do begin
    if i ne except then begin
      this_inside=(*(self.elements))[i]->is_inside(r);
      if self.is_union then begin
        result=logical_or(result,this_inside)
      end else begin
        result=logical_and(result,this_inside)
      end
    end
  end
  return, result
end

;function intersect
;  Calculates intersection parameter t for this surface and a ray
;  R=R0+V*t
;Inputs:
;  R0 - Origin of ray to test. Must be a single vector
;  V  - Direction of ray to test. May be a grid
;Returns:
;  t parameter of intersection, or negative value or NaN if no intersection.
function csg::intersect,r0,v
  ;The basic idea is to find the intersection with all the elements and choose
  ;the "best" valid intersection. The best one is closest.

  ;Start with a "best" intersection which is very far away
  t=dblarr(shape_grid(v))+9d99

  ;Test this ray against all the elements
  for i=0,n_elements((*(self.elements)))-1 do begin
    ;Calculate the intersection with this element
    this_intersect=(*(self.elements))[i]->intersect(r0,v)

    ;We only care about rays where this intersect is both better and valid
    is_intersect=logical_and(this_intersect lt t,this_intersect ge 0)
    w_is_intersect=where(is_intersect,num_is_intersect)
    if num_is_intersect gt 0 then begin
      ;For those rays are better and valid, check if it satisfies the set relation
      if self.is_union then begin
        ;For a union, it's easy. If it intersects, it counts
        t[w_is_intersect]=this_intersect[w_is_intersect]
      end else begin
        ;For an intersection, it only counts if the point is inside all the other
        ;elements. We already know it's exactly on the surface of this element
        is_i2sect=logical_and(is_intersect,self->is_inside_except(r0+v*this_intersect,i))
        w_is_i2sec=where(is_i2sect,num_is_i2sect)
        if num_is_i2sect gt 0 then begin
          t[w_is_i2sect]=this_intersect[w_is_i2sect]
        end
      end
    end
  end
  ;Clean up those places where there never was an intersect and the "best" is still 9d99
  w_no_intersect=where(t gt 1d99,n_no_intersect)
  if n_no_intersect gt 0 then t[w_no_intersect]=-1

  return,t
end

; This method will fail unless the csg elements have their own 'to_pov' function
function csg::to_pov
  if self.is_union then begin
    result="union {"+string([13b,10b])
  end else begin
    result="intersection {"+string([13b,10b])
  end
  for i=0,n_elements((*(self.elements)))-1 do begin
    result=result+(*(self.elements))[i]->to_pov()+string([13b,10b])
  end
  result=result+"}"
  return,result
end
function csg::to_pov_array
  result="{"+string([13b,10b])
  for i=0,n_elements((*(self.elements)))-1 do begin
    result=result+(*(self.elements))[i]->to_pov_array()+string([13b,10b])
  end
  result=result+"},"+string([13b,10b])
  return,result
end

;--------------------------------------------------------------------------
; CLEANUP
; --------------------------------------------------------------------------
; PURPOSE:
;   Lifecycle method to clean itself up.
;
pro csg::CLEANUP
  str_free,self.elements
  self->surface::cleanup
end

; --------------------------------------------------------------------------
; Object Data Structure Definition
; --------------------------------------------------------------------------
; PURPOSE:
;   Define the data structure associated with this object.

pro csg__define
  self = {csg, inherits surface,is_union:0,elements:ptr_new()}
end