Ruby: Construction of a closed curve with discrete edge points, that calculate its area and the center of mass.

gistfile1.rb

class Point
  def initialize(x, y)
    self.x = x
    self.y = y
  end

  def inspect
    "(#{x}, #{y})"
  end

  attr_accessor :x, :y
end


# Do not treat a line perpendicular to x axis
class LineSegment

  attr_accessor :from, :to
  attr_accessor :from_x, :to_x # X dimension order
  attr_accessor :from_y, :to_y # Y dimension order
  attr_accessor :angle, :angle_inverse

  def initialize(one, another)
    set_points(one, another)
    set_angle
  end

  def set_points(one, another)
    if one.x == another.x
      one.x+=1e-10
    end
    self.from , self.to = one, another

    self.from_x, self.to_x = [one, another].sort{|a,b|a.x<=>b.x}
    self.from_y, self.to_y = [one, another].sort{|a,b|a.y<=>b.y}
  end

  def set_angle
    self.angle = (self.to.y - self.from.y) / (self.to.x - self.from.x)
    self.angle_inverse = 1.0 / self.angle
  end

  def y_for(x)
    self.from_y.y + self.angle * (x - self.from_y.x)
  end

  def x_for(y)
    self.from_x.x + self.angle_inverse * (y - self.from_x.y)
  end

  def includes_x(x)
    (self.from_x.x <= x && x <= self.to_x.x)
  end

  def includes_y(y)
    (self.from_y.y <= y && y <= self.to_y.y)
  end

  def inspect
    "#{self.from} <-> #{self.to}"
  end

end

class ApproximatedClosedCurve
  attr_accessor :line_segments
  attr_accessor :calculation_mesh
  attr_accessor :bottom_left_point, :top_right_point
  attr_accessor :points
  def initialize
    self.line_segments = []
    self.calculation_mesh = 1.0
    area = 0
    area_calculated = false
  end

  # Requires at least 3 points
  def construct(points)
    self.points = points
    area_calculated = false
    points.each_with_index do |point, i|
      next_point = i==points.size-1 ? points[0] : points[i+1]
      line_segment = LineSegment.new(point, next_point)
      line_segments << line_segment
    end

    calculate_range
  end

  def calculate_range
    min_x = points.min { |a, b| a.x <=> b.x}.x
    min_y = points.min { |a, b| a.y <=> b.y}.y
    max_x = points.max { |a, b| a.x <=> b.x}.x
    max_y = points.max { |a, b| a.y <=> b.y}.y
    self.bottom_left_point = Point.new(min_x, min_y)
    self.top_right_point = Point.new(max_x, max_y)
  end

  def area
    weighted_integral_over_x_direction(->(x){1.0})
  end

  def weighted_integral_over_x_direction(weight_function)
    current_x = bottom_left_point.x + 0.0001
    res_x = 0
    while current_x + calculation_mesh < top_right_point.x
      #STDOUT.puts "current X #{current_x} : area: #{res}"
      next_x = current_x + calculation_mesh
      one_line_segments = line_segments_include_x(current_x)
      another_line_segments = line_segments_include_x(next_x)
      pairs = [one_line_segments.size, another_line_segments.size].min / 2
      i = 0
      pairs.times do |pair|
        left_pair = [one_line_segments[i], one_line_segments[i+1]]
        right_pair = [another_line_segments[i], another_line_segments[i+1]]
        left_y_diff = (left_pair[0].y_for(current_x) - left_pair[1].y_for(current_x)).abs
        right_y_diff = (right_pair[0].y_for(next_x) - right_pair[1].y_for(next_x)).abs
        this_val = calculation_mesh * (left_y_diff + right_y_diff) / 2.0 * weight_function.call(current_x)
        res_x += this_val
        i+=2
      end
      current_x += calculation_mesh
    end
    res_x
  end

  def weighted_integral_over_y_direction(weight_function)
    current_y = bottom_left_point.y + 0.0001
    res_y = 0
    while current_y + calculation_mesh < top_right_point.y
#      STDOUT.puts "current Y #{current_y} : area: #{res_y}"
      next_y = current_y + calculation_mesh
      one_line_segments = line_segments_include_y(current_y)
      another_line_segments = line_segments_include_y(next_y)
      pairs = [one_line_segments.size, another_line_segments.size].min / 2
      i = 0
      pairs.times do |pair|
        lower_pair = [one_line_segments[i], one_line_segments[i+1]]
        upper_pair = [another_line_segments[i], another_line_segments[i+1]]
        lower_x_diff = (lower_pair[0].x_for(current_y) - lower_pair[1].x_for(current_y)).abs
        upper_x_diff = (upper_pair[0].x_for(next_y) - upper_pair[1].x_for(next_y)).abs
        this_val = calculation_mesh * (lower_x_diff + upper_x_diff) / 2.0 * weight_function.call(current_y)
        res_y += this_val
        i+=2
      end
      current_y += calculation_mesh
    end
    res_y
  end

  def mass_center
    res_x = weighted_integral_over_x_direction(->(x){x + calculation_mesh/2.0})
    res_y = weighted_integral_over_y_direction(->(y){y + calculation_mesh/2.0})
    _area = self.area
    Point.new(res_x/_area, res_y/_area)
  end

  def line_segments_include_x(x)
    line_segments.select{|s|s.includes_x(x)}
  end

  def line_segments_include_y(y)
    line_segments.select{|s|s.includes_y(y)}
  end

end