Bisect a polygon by angle
arcgissamples\geometry\BisectPolygonByAngle.java
/* Copyright 2010 ESRI
* 
* All rights reserved under the copyright laws of the United States
* and applicable international laws, treaties, and conventions.
* 
* You may freely redistribute and use this sample code, with or
* without modification, provided you include the original copyright
* notice and use restrictions.
* 
* See the use restrictions.
* 
*/
package arcgissamples.geometry;

import java.io.*;
import java.net.UnknownHostException;

import com.esri.arcgis.datasourcesfile.*;
import com.esri.arcgis.geodatabase.*;
import com.esri.arcgis.geometry.*;
import com.esri.arcgis.interop.AutomationException;
import com.esri.arcgis.system.*;

public class BisectPolygonByAngle
{
  public GeometryBag BisectByAngle(Polygon inputPolygon, double bisectRatio, double inputAngle, boolean isLeft,
      double precision) throws Exception
  {

    boolean isHorizontal = false;
    double ratio = 0.0, leftPolygonArea = 0.0;
    int loopCtr = 0;
    IGeometry[] leftGeometry = new IGeometry[1];
    IGeometry[] rightGeometry = new IGeometry[1];
    double pi = Math.PI;

    if (inputPolygon == null)
    {
      throw new Exception("The polygon is not valid");
    }
    if (inputPolygon.isEmpty())
    {
      throw new Exception("The polygon is not empty");
    }
    // Check if the angle is between 0 and 2PI
    if (!(0 <= inputAngle) || !(inputAngle <= (2 * pi)))
    {
      throw new Exception("Angle should be between 0 and 2PI Radian");
    }
    // Check if the angle is 0
    if (inputAngle == 0 || inputAngle == 2 * pi || inputAngle == -2 * pi)
    {
      isHorizontal = true;
    }
    // If input angle is larger than PI then take the corresponding angle in opposite quadrant
    if (inputAngle > pi)
    {
      inputAngle -= pi;
    }

    // Check if the polygon as less than two exterior rings, cannot handle more than that
    if (inputPolygon.getExteriorRingCount() > 1)
    {
      throw new Exception("Select a polygon with less than two exterior rings.");
    }

    // Get polygon envelope
    Envelope inputPolygonEnvelope = (Envelope) inputPolygon.getEnvelope();

    // Get the Area of the polygon
    double fullArea = inputPolygon.getArea();

    // Create the baselines, the baselines are used as start points in the iterative process.
    GeometryBag baselines = createBaseLineByAngle(inputPolygonEnvelope, inputAngle, isHorizontal);

    // Get the min and max lines
    Line minLine = (Line) baselines.getGeometry(0);
    Line maxLine = (Line) baselines.getGeometry(1);

    int maxIterations = 1000; // That value could be modified
    // Reverse the ratio if not "left"
    if (!isLeft)
    {
      bisectRatio = 1 - bisectRatio;
    }
    while (ratio != bisectRatio)
    {
      loopCtr++;
      // Safety check to avoid infinite loop
      if (loopCtr == maxIterations)
      {
        throw new Exception("Cannot determine the bisection !");
      }
      // Create the line in the middle of the min and max lines
      Line halfLine = createHalfLineByAngle(minLine, maxLine);

      // Wrap the created line into a polyline to cut the polygon
      Polyline polylineCutter = new Polyline();
      polylineCutter.addSegment(halfLine, null, null);

      // The cutter will be clipped to the extent of the polygon's
      // extent
      Envelope polyEnvelope = (Envelope) inputPolygon.getEnvelope();
      polyEnvelope.expand(1.01, 1.01, true);

      // Clip the cutter
      polylineCutter.clip(polyEnvelope);

      // Make sure the cutter crosses the polylgon
      if (inputPolygon.crosses(polylineCutter))
      {

        // Cut the polygon
        inputPolygon.cut(polylineCutter, leftGeometry, rightGeometry);
        if ((leftGeometry == null) || (rightGeometry == null))
        {
          return null;
        }
        if (leftGeometry[0].isEmpty() || rightGeometry[0].isEmpty())
        {
          return null;
        }
        Polygon leftPolygon = (Polygon) leftGeometry[0];
        leftPolygonArea = leftPolygon.getArea();
        ratio = leftPolygonArea / fullArea;

      }
      else
      {
        // When the cutter doesn't cross the polygon, define the ratio
        if (bisectRatio > 0.5)
        {
          // To make sure the min and max line are determined correctly
          ratio = 1.0;
        }
        else
        {
          ratio = 0.0;
        }
      }
      // Check if the cut polygons meet the precision
      if (Math.abs(ratio - bisectRatio) < precision)
      {
        break;
      }
      if (ratio < bisectRatio)
      {
        minLine = halfLine;
      }
      else
      {
        maxLine = halfLine;
      }
    }
    // Create the output geometry bag
    GeometryBag geometryBag = new GeometryBag();
    geometryBag.setSpatialReferenceByRef(inputPolygon.getSpatialReference());

    geometryBag.addGeometry(leftGeometry[0], null, null);
    geometryBag.addGeometry(rightGeometry[0], null, null);

    Polygon leftPolygon = (Polygon) leftGeometry[0];
    Polygon rightPolygon = (Polygon) rightGeometry[0];

    if (!withinPrecision(leftPolygon.getArea() + rightPolygon.getArea(), fullArea, precision))
    {
      throw new Exception("Cannot precisely determine the bisection. Try reducing the precision.");
    }

    return geometryBag;
  }

  private Line createHalfLineByAngle(ILine pMinLine, ILine pMaxLine) throws Exception
  {
    Line line = new Line();
    Point point0 = new Point();
    Point point1 = new Point();
    point0.putCoords((pMinLine.getFromPoint().getX() + pMaxLine.getFromPoint().getX()) / 2, (pMinLine
        .getFromPoint().getY() + pMaxLine.getFromPoint().getY()) / 2);
    point1.putCoords((pMinLine.getToPoint().getX() + pMaxLine.getToPoint().getX()) / 2, (pMinLine.getToPoint()
        .getY() + pMaxLine.getToPoint().getY()) / 2);
    line.putCoords(point0, point1);
    return line;
  }

  private GeometryBag createBaseLineByAngle(Envelope envelope, double angle, boolean isHorizontal) throws Exception
  {
    double pi = Math.PI;
    Line minLine, maxLine, lftBaseline;
    Point maxPoint1, maxPoint2, minPoint1, minPoint2, lftpb0, lftpb1, midPoint0, midPoint1;
    Point constructPointMin, constructPointMax, constructPointMid0, constructPointMid1;
    GeometryBag geometryBag;
    double x0, y0, x1, y1;
    ITransform2D transform2D0, transform2D1;
    if (isHorizontal)
    {
      // If the input angle is 0 the maxline and min lines
      // are equals to the edges of the envelope
      minLine = new Line();
      maxLine = new Line();
      maxPoint1 = new Point();
      maxPoint2 = new Point();
      minPoint1 = new Point();
      minPoint2 = new Point();
      maxPoint2.putCoords(envelope.getXMin(), envelope.getYMax());
      maxPoint1.putCoords(envelope.getXMax(), envelope.getYMax());
      minPoint2.putCoords(envelope.getXMin(), envelope.getYMin());
      minPoint1.putCoords(envelope.getXMax(), envelope.getYMin());
      maxLine.putCoords(maxPoint1, maxPoint2);
      minLine.putCoords(minPoint1, minPoint2);
    }
    else
    {
      // If the angle is not horizontal then the
      // 'baselines are based on the egdes of the envelope +/- the input angle
      minLine = new Line();
      maxLine = new Line();
      maxPoint1 = new Point();
      minPoint2 = new Point();
      minPoint2.putCoords(envelope.getXMin(), envelope.getYMax());
      constructPointMin = new Point();
      lftBaseline = new Line();
      lftpb0 = new Point();
      lftpb0.putCoords(envelope.getXMin(), envelope.getYMax());
      lftpb1 = new Point();
      lftpb1.putCoords(envelope.getXMin(), envelope.getYMin());
      lftBaseline.putCoords(lftpb0, lftpb1);
      constructPointMin.constructDeflection(lftBaseline, 2 * (lftBaseline.getLength() / Math
          .cos((pi / 2) - angle)), -((pi / 2) - angle));
      minPoint1 = constructPointMin;
      minLine.putCoords(minPoint1, minPoint2);
      maxPoint1.putCoords(envelope.getXMax(), envelope.getYMin());
      constructPointMax = new Point();
      constructPointMax.constructAngleDistance(maxPoint1, angle, minLine.getLength());
      maxPoint2 = constructPointMax;
      maxLine.putCoords(maxPoint1, maxPoint2);
    }
    geometryBag = new GeometryBag();
    if (angle > pi / 2)
    {
      // In this case have to move the line
      constructPointMid0 = new Point();
      constructPointMid0.constructAlong(minLine, esriSegmentExtension.esriNoExtension, 0.5, true);
      midPoint0 = constructPointMid0;
      x0 = envelope.getXMin() - midPoint0.getX();
      y0 = envelope.getYMin() - midPoint0.getY();
      transform2D0 = minLine;
      transform2D0.move(x0, y0);
      geometryBag.addGeometry((IGeometry) transform2D0, null, null);
      constructPointMid1 = new Point();
      constructPointMid1.constructAlong(maxLine, esriSegmentExtension.esriNoExtension, 0.5, true);
      midPoint1 = constructPointMid1;
      x1 = envelope.getXMax() - midPoint1.getX();
      y1 = envelope.getYMax() - midPoint1.getY();
      transform2D1 = maxLine;
      transform2D1.move(x1, y1);
      geometryBag.addGeometry((IGeometry) transform2D1, null, null);
    }
    else
    {
      geometryBag.addGeometry(minLine, null, null);
      geometryBag.addGeometry(maxLine, null, null);
    }
    return geometryBag;
  }

  private boolean withinPrecision(double TestValue, double ExpectedValue, double Precision)
  {
    double preciseTestValue = Math.floor(TestValue * (1.0 / Precision));
    double preciseExpectedValue = Math.floor(ExpectedValue * (1.0 / Precision));
    if (Math.abs(preciseTestValue - preciseExpectedValue) > 1)
    {
      return false;
    }
    return true;
  }

  public static void main(String[] args)
  {
    if (args.length != 1)
    {
      System.out.println("Usage: BisetAPolygonByAngle <Source ShapeFile>");
      System.exit(0);
    }

    String inShape = args[0];
    File shapefileFile = new File(inShape);
    if (!shapefileFile.exists())
    {
      System.out.println("Doesn't exist: " + shapefileFile.getAbsolutePath());
      System.exit(0);
    }

    EngineInitializer.initializeEngine();
    
    AoInitialize aoInit = null;
    try
    {
      aoInit = new AoInitialize();
      initializeArcGISLicenses(aoInit);
    }
    catch (UnknownHostException e1)
    {
      e1.printStackTrace();
    }
    catch (IOException e1)
    {
      e1.printStackTrace();
    }
    
    // Open a shapefile to get a polygon...
    GeometryBag result = null;
    try
    {
      ShapefileWorkspaceFactory shapefileWorkspaceFactory = new ShapefileWorkspaceFactory();
      Workspace iWorkspace = (Workspace) shapefileWorkspaceFactory.openFromFile(shapefileFile.getParent(), 0);
      IFeatureClass iFeatureClass = iWorkspace.openFeatureClass(shapefileFile.getName());

      // You may need to change the following index, if that polygon has more than 1 ring...
      IFeature iFeature = iFeatureClass.getFeature(1);

      // If the geometry type is not of type Polygon, exit...
      if (iFeature.getShape().getGeometryType() != esriGeometryType.esriGeometryPolygon)
      {
        System.out.println("The input shapefile must be of type Polygon.");
        try
        {
          aoInit.shutdown();
        }
        catch (Exception e)
        {
          e.printStackTrace();
        }
        System.exit(0);
      }

      Polygon poly = (Polygon) iFeature.getShape();

      System.out.println("");
      System.out.println("Area of input Polygon: " + poly.getArea());

      // Call the method that performs the bisection...
      BisectPolygonByAngle bisect = new BisectPolygonByAngle();
      result = bisect.BisectByAngle(poly, 0.0, 0.9, false, 0.01);

      // Show the results...
      System.out.println("");
      int numPolyResult = result.getGeometryCount();
      System.out.println("Number of resulting Polygons from the split: " + numPolyResult);

      for (int i = 0; i < numPolyResult; i++)
      {
        IArea area = (IArea) result.getGeometry(i);
        System.out.println("Area of Polygon " + (i + 1) + ": " + area.getArea());
      }

      System.out.println();
      System.out.println("Done!");
      
      aoInit.shutdown();

    }
    catch (Exception ex)
    {
      ex.printStackTrace();
      try
      {
        aoInit.shutdown();
      }
      catch (AutomationException e)
      {
        e.printStackTrace();
      }
      catch (IOException e)
      {
        e.printStackTrace();
      }
      System.exit(0);
    }
  }
  
  private static void initializeArcGISLicenses(AoInitialize aoInit) {
    try {
      if (aoInit.isProductCodeAvailable(esriLicenseProductCode.esriLicenseProductCodeEngine) 
          == esriLicenseStatus.esriLicenseAvailable)
        aoInit.initialize(esriLicenseProductCode.esriLicenseProductCodeEngine);
      else if (aoInit.isProductCodeAvailable(esriLicenseProductCode.esriLicenseProductCodeArcView) 
          == esriLicenseStatus.esriLicenseAvailable)
        aoInit.initialize(esriLicenseProductCode.esriLicenseProductCodeArcView);
      else if (aoInit.isProductCodeAvailable(esriLicenseProductCode.esriLicenseProductCodeArcEditor) 
          == esriLicenseStatus.esriLicenseAvailable)
        aoInit.initialize(esriLicenseProductCode.esriLicenseProductCodeArcEditor);
      else if (aoInit.isProductCodeAvailable(esriLicenseProductCode.esriLicenseProductCodeArcInfo) == com.esri.arcgis.system.esriLicenseStatus.esriLicenseAvailable)
        aoInit.initialize(esriLicenseProductCode.esriLicenseProductCodeArcInfo);
      else{
        System.err.println("Could not initialize an Engine or ArcView license. Exiting application.");
        System.exit(-1);
      }
    } catch (Exception e) {e.printStackTrace();}
  }

}