overlayop.cpp

来自「一个很好的vc底层代码」· C++ 代码 · 共 941 行 · 第 1/2 页

 * They do not form part of the result area boundary.
 */
void
OverlayOp::findResultAreaEdges(int opCode)
{
	vector<EdgeEnd*> *ee=graph->getEdgeEnds();
	for(unsigned int i=0;i<ee->size();i++) {
		DirectedEdge *de=(DirectedEdge*) (*ee)[i];
		// mark all dirEdges with the appropriate label
		Label *label=de->getLabel();
		if (label->isArea()
			&& !de->isInteriorAreaEdge()
			&& isResultOfOp(label->getLocation(0,Position::RIGHT),
							label->getLocation(1,Position::RIGHT),
							opCode)
			) {
				de->setInResult(true);
				//Debug.print("in result "); Debug.println(de);
		}
	}
}

/*
 * If both a dirEdge and its sym are marked as being in the result, cancel
 * them out.
 */
void
OverlayOp::cancelDuplicateResultEdges()
{
	// remove any dirEdges whose sym is also included
	// (they "cancel each other out")
	vector<EdgeEnd*> *ee=graph->getEdgeEnds();
	for(int i=0;i<(int)ee->size();i++) {
		DirectedEdge *de=(DirectedEdge*) (*ee)[i];
		DirectedEdge *sym=de->getSym();
		if (de->isInResult() && sym->isInResult()) {
			de->setInResult(false);
			sym->setInResult(false);
			//Debug.print("cancelled "); Debug.println(de); Debug.println(sym);
		}
	}
}

/*
 * This method is used to decide if a point node should be included
 * in the result or not.
 *
 * @return true if the coord point is covered by a result Line or
 *	Area geometry
 */
bool
OverlayOp::isCoveredByLA(const Coordinate& coord)
{
	if (isCovered(coord,resultLineList)) return true;
	if (isCovered(coord,resultPolyList)) return true;
	return false;
}

/*
 * This method is used to decide if an L edge should be included
 * in the result or not.
 *
 * @return true if the coord point is covered by a result Area geometry
 */
bool
OverlayOp::isCoveredByA(const Coordinate& coord)
{
	if (isCovered(coord,resultPolyList)) return true;
	return false;
}

/*
 * @return true if the coord is located in the interior or boundary of
 * a geometry in the list.
 */
bool
OverlayOp::isCovered(const Coordinate& coord,vector<Geometry*> *geomList)
{
	for(int i=0;i<(int)geomList->size();i++) {
		Geometry *geom=(*geomList)[i];
		int loc=ptLocator->locate(coord,geom);
		if (loc!=Location::EXTERIOR) return true;
	}
	return false;
}

/*
 * @return true if the coord is located in the interior or boundary of
 * a geometry in the list.
 */
bool
OverlayOp::isCovered(const Coordinate& coord,vector<LineString*> *geomList)
{
	for(int i=0;i<(int)geomList->size();i++) {
		Geometry *geom=(Geometry*)(*geomList)[i];
		int loc=ptLocator->locate(coord,geom);
		if (loc!=Location::EXTERIOR) return true;
	}
	return false;
}

/*
 * @return true if the coord is located in the interior or boundary of
 * a geometry in the list.
 */
bool
OverlayOp::isCovered(const Coordinate& coord,vector<Polygon*> *geomList)
{
	for(int i=0;i<(int)geomList->size();i++) {
		Geometry *geom=(Geometry*)(*geomList)[i];
		int loc=ptLocator->locate(coord,geom);
		if (loc!=Location::EXTERIOR) return true;
	}
	return false;
}

/*
 * Build a Geometry containing given elements
 * This function will take ownership of vectors elements.
 */
Geometry*
OverlayOp::computeGeometry(vector<Point*> *nResultPointList,
                              vector<LineString*> *nResultLineList,
                              vector<Polygon*> *nResultPolyList) {
	int i;
	vector<Geometry*> *geomList=new vector<Geometry*>();
	// element geometries of the result are always in the order P,L,A
	for(i=0;i<(int)nResultPointList->size();i++) {
		Point *pt=(*nResultPointList)[i];
		geomList->push_back(pt);
	}
	for(i=0;i<(int)nResultLineList->size();i++) {
		LineString *ls=(*nResultLineList)[i];
		geomList->push_back(ls);
	}
	for(i=0;i<(int)nResultPolyList->size();i++) {
		Polygon *q=(*nResultPolyList)[i];
		geomList->push_back(q);
	}
	// build the most specific geometry possible
	Geometry *g=geomFact->buildGeometry(geomList);
	return g;
}

void
OverlayOp::computeOverlay(int opCode)
	//throw(TopologyException *)
{
	vector<Edge*> *baseSplitEdges=NULL;
	PolygonBuilder *polyBuilder=NULL;
	LineBuilder *lineBuilder=NULL;
	PointBuilder *pointBuilder=NULL;

	// copy points from input Geometries.
	// This ensures that any Point geometries
	// in the input are considered for inclusion in the result set

	copyPoints(0);
	copyPoints(1);


	// node the input Geometries
	delete (*arg)[0]->computeSelfNodes(li,false);
	delete (*arg)[1]->computeSelfNodes(li,false);

#if DEBUG
	cerr<<"OverlayOp::computeOverlay: computed SelfNodes"<<endl;
#endif

	// compute intersections between edges of the two input geometries
	delete (*arg)[0]->computeEdgeIntersections((*arg)[1],li,true);

#if DEBUG
	cerr<<"OverlayOp::computeOverlay: computed EdgeIntersections"<<endl;
	cerr<<"OverlayOp::computeOverlay: li: "<<li->toString()<<endl;
#endif


	baseSplitEdges = new vector<Edge*>();
	(*arg)[0]->computeSplitEdges(baseSplitEdges);
	(*arg)[1]->computeSplitEdges(baseSplitEdges);

	// add the noded edges to this result graph
	insertUniqueEdges(baseSplitEdges);
	computeLabelsFromDepths();
	replaceCollapsedEdges();
	//Debug.println(edgeList);

    	// debugging only
    	//NodingValidator nv = new NodingValidator(edgeList.getEdges());
    	//nv.checkValid();

	graph->addEdges(edgeList->getEdges());

	try {
		// this can throw TopologyException *
		computeLabelling();

		//Debug.printWatch();
		labelIncompleteNodes();
		//Debug.printWatch();
		//nodeMap.print(System.out);


		/*
		 * The ordering of building the result Geometries is important.
		 * Areas must be built before lines, which must be built
		 * before points.
		 * This is so that lines which are covered by areas are not
		 * included explicitly, and similarly for points.
		 */
		findResultAreaEdges(opCode);
		cancelDuplicateResultEdges();

		polyBuilder=new PolygonBuilder(geomFact,cga);
		
		// might throw a TopologyException *
		polyBuilder->add(graph);

		vector<Geometry*> *gv=polyBuilder->getPolygons();
		resultPolyList=new vector<Polygon*>();
		for(int i=0;i<(int)gv->size();i++) {
			resultPolyList->push_back((Polygon*)(*gv)[i]);
		}
		delete gv;
		lineBuilder=new LineBuilder(this,geomFact,ptLocator);
		resultLineList=lineBuilder->build(opCode);
		pointBuilder=new PointBuilder(this,geomFact,ptLocator);
		resultPointList=pointBuilder->build(opCode);
		// gather the results from all calculations into a single
		// Geometry for the result set
		resultGeom=computeGeometry(resultPointList,resultLineList,resultPolyList);
#if USE_ELEVATION_MATRIX
		elevationMatrix->elevate(resultGeom);
#endif // USE_ELEVATION_MATRIX
		


	} catch (...) {
		delete baseSplitEdges;
		delete polyBuilder;
		delete lineBuilder;
		delete pointBuilder;
		throw;
	}


	delete polyBuilder;
	delete lineBuilder;
	delete pointBuilder;
	delete baseSplitEdges;
}

/*
 * Insert an edge from one of the noded input graphs.
 * Checks edges that are inserted to see if an
 * identical edge already exists.
 * If so, the edge is not inserted, but its label is merged
 * with the existing edge.
 */
void
OverlayOp::insertUniqueEdge(Edge *e)
{
	//Debug.println(e);
#if DEBUG
	cerr<<"OverlayOp::insertUniqueEdge("<<e->print()<<")"<<endl;
#endif

	int foundIndex=edgeList->findEdgeIndex(e);
	// If an identical edge already exists, simply update its label
	if (foundIndex>=0) {
#if DEBUG
		cerr<<"  found identical edge, should merge Z"<<endl;
#endif
		Edge *existingEdge=edgeList->get(foundIndex);
		Label *existingLabel=existingEdge->getLabel();
		Label *labelToMerge=e->getLabel();

		// check if new edge is in reverse direction to existing edge
		// if so, must flip the label before merging it
		if (!existingEdge->isPointwiseEqual(e)) {
//			labelToMerge=new Label(e->getLabel());
			labelToMerge->flip();
		}
		Depth *depth=existingEdge->getDepth();
		// if this is the first duplicate found for this edge, initialize the depths
		///*
		if (depth->isNull()) {
			depth->add(existingLabel);
		}
		//*/
		depth->add(labelToMerge);

		existingLabel->merge(labelToMerge);
		//Debug.print("inserted edge: "); Debug.println(e);
		//Debug.print("existing edge: "); Debug.println(existingEdge);
		dupEdges.push_back(e);
	} else {  // no matching existing edge was found
#if DEBUG
		cerr<<"  no matching existing edge"<<endl;
#endif
		// add this new edge to the list of edges in this graph
		//e.setName(name+edges.size());
		//e.getDepth().add(e.getLabel());
		edgeList->add(e);
	}
}

/*
 * Update the labels for edges according to their depths.
 * For each edge, the depths are first normalized.
 * Then, if the depths for the edge are equal,
 * this edge must have collapsed into a line edge.
 * If the depths are not equal, update the label
 * with the locations corresponding to the depths
 * (i.e. a depth of 0 corresponds to a Location of EXTERIOR,
 * a depth of 1 corresponds to INTERIOR)
 */
void
OverlayOp::computeLabelsFromDepths()
{
	for(int j=0;j<(int)edgeList->getEdges()->size();j++) {
		Edge *e=edgeList->get(j);
		Label *lbl=e->getLabel();
		Depth *depth=e->getDepth();
		/**
		 * Only check edges for which there were duplicates,
		 * since these are the only ones which might
		 * be the result of dimensional collapses.
		 */
		if (!depth->isNull()) {
			depth->normalize();
			for (int i=0;i<2;i++) {
				if (!lbl->isNull(i) && lbl->isArea() && !depth->isNull(i)) {
					/**
					* if the depths are equal, this edge is the result of
					* the dimensional collapse of two or more edges.
					* It has the same location on both sides of the edge,
					* so it has collapsed to a line.
					*/
					if (depth->getDelta(i)==0) {
						lbl->toLine(i);
					} else {
						/**
						* This edge may be the result of a dimensional collapse,
						* but it still has different locations on both sides.  The
						* label of the edge must be updated to reflect the resultant
						* side locations indicated by the depth values.
						*/
						Assert::isTrue(!depth->isNull(i,Position::LEFT),"depth of LEFT side has not been initialized");
						lbl->setLocation(i,Position::LEFT,depth->getLocation(i,Position::LEFT));
						Assert::isTrue(!depth->isNull(i,Position::RIGHT),"depth of RIGHT side has not been initialized");
						lbl->setLocation(i,Position::RIGHT,depth->getLocation(i,Position::RIGHT));
					}
				}
			}
		}
	}
}

} // namespace geos

/**********************************************************************
 * $Log: OverlayOp.cpp,v $
 * Revision 1.37.2.1  2005/05/23 16:41:45  strk
 * Back-ported removal of useless Coordinate copies in mergeZ() - patch by Safe Software
 *
 * Revision 1.37  2004/12/08 14:31:17  strk
 * elevationMatrix deleted by destructor
 *
 * Revision 1.36  2004/12/08 13:54:44  strk
 * gcc warnings checked and fixed, general cleanups.
 *
 * Revision 1.35  2004/11/29 16:05:33  strk
 * Fixed a bug in LineIntersector::interpolateZ causing NaN values
 * to come out.
 * Handled dimensional collapses in ElevationMatrix.
 * Added ISNAN macro and changed ISNAN/FINITE macros to avoid
 * dispendious isnan() and finite() calls.
 *
 * Revision 1.34  2004/11/26 09:22:50  strk
 * Added FINITE(x) macro and its use.
 * Made input geoms average Z computation optional in OverlayOp.
 *
 * Revision 1.33  2004/11/24 18:10:42  strk
 * Stricter handling of USE_ELEVATION_MATRIX define
 *
 * Revision 1.32  2004/11/23 19:53:07  strk
 * Had LineIntersector compute Z by interpolation.
 *
 * Revision 1.31  2004/11/23 16:22:49  strk
 * Added ElevationMatrix class and components to do post-processing draping of overlayed geometries.
 *
 * Revision 1.30  2004/11/22 15:51:52  strk
 * Added interpolation of containing geometry's average Z for point_in_poly case.
 *
 * Revision 1.29  2004/11/22 11:34:49  strk
 * More debugging lines and comments/indentation cleanups
 *
 * Revision 1.28  2004/11/20 17:16:10  strk
 * Handled Z merging for point on polygon boundary case.
 *
 * Revision 1.27  2004/11/20 16:25:17  strk
 * Added Z computation for point on line case.
 *
 * Revision 1.26  2004/11/17 15:09:08  strk
 * Changed COMPUTE_Z defaults to be more conservative
 *
 * Revision 1.25  2004/11/17 08:13:16  strk
 * Indentation changes.
 * Some Z_COMPUTATION activated by default.
 *
 * Revision 1.24  2004/10/21 22:29:54  strk
 * Indentation changes and some more COMPUTE_Z rules
 *
 * Revision 1.23  2004/10/20 17:32:14  strk
 * Initial approach to 2.5d intersection()
 *
 * Revision 1.22  2004/07/02 13:28:29  strk
 * Fixed all #include lines to reflect headers layout change.
 * Added client application build tips in README.
 *
 * Revision 1.21  2004/07/01 14:12:44  strk
 *
 * Geometry constructors come now in two flavors:
 * 	- deep-copy args (pass-by-reference)
 * 	- take-ownership of args (pass-by-pointer)
 * Same functionality is available through GeometryFactory,
 * including buildGeometry().
 *
 * Revision 1.20  2004/06/15 20:13:42  strk
 * updated to respect deep-copy GeometryCollection interface
 *
 * Revision 1.19  2004/05/17 08:34:31  strk
 * reduced stack allocations, try/catch blocks in ::overlayOp
 *
 * Revision 1.18  2004/05/03 10:43:43  strk
 * Exception specification considered harmful - left as comment.
 *
 * Revision 1.17  2004/04/20 13:24:15  strk
 * More leaks removed.
 *
 * Revision 1.16  2004/04/14 13:14:29  strk
 * Removed deletion of externally pointed GeometryFactory from OverlayOp destructor
 *
 * Revision 1.15  2004/04/10 08:40:01  ybychkov
 * "operation/buffer" upgraded to JTS 1.4
 *
 * Revision 1.14  2004/03/29 06:59:25  ybychkov
 * "noding/snapround" package ported (JTS 1.4);
 * "operation", "operation/valid", "operation/relate" and "operation/overlay" upgraded to JTS 1.4;
 * "geom" partially upgraded.
 *
 * Revision 1.13  2004/03/19 09:48:46  ybychkov
 * "geomgraph" and "geomgraph/indexl" upgraded to JTS 1.4
 *
 * Revision 1.12  2003/11/12 18:02:56  strk
 * Added throw specification. Fixed leaks on exceptions.
 *
 * Revision 1.11  2003/11/12 16:14:56  strk
 * Added some more throw specifications and cleanup on exception (leaks removed).
 *
 * Revision 1.10  2003/11/07 01:23:42  pramsey
 * Add standard CVS headers licence notices and copyrights to all cpp and h
 * files.
 *
 *
 **********************************************************************/