Warning, file /education/marble/tools/vectorosm-tilecreator/VectorClipper.cpp was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).
0001 // SPDX-License-Identifier: LGPL-2.1-or-later 0002 // 0003 // SPDX-FileCopyrightText: 2016 Dennis Nienhüser <nienhueser@kde.org> 0004 // SPDX-FileCopyrightText: 2016 David Kolozsvari <freedawson@gmail.com> 0005 // 0006 0007 #include "VectorClipper.h" 0008 0009 #include "TileId.h" 0010 0011 #include "GeoDataLatLonAltBox.h" 0012 #include "GeoDataPolygon.h" 0013 #include "GeoDataRelation.h" 0014 #include "OsmObjectManager.h" 0015 #include "TileCoordsPyramid.h" 0016 0017 0018 #include <QDebug> 0019 #include <QPolygonF> 0020 #include <QPair> 0021 #include <QStringBuilder> 0022 0023 namespace Marble { 0024 0025 VectorClipper::VectorClipper(GeoDataDocument* document, int maxZoomLevel) : 0026 m_maxZoomLevel(maxZoomLevel) 0027 { 0028 for (auto feature: document->featureList()) { 0029 if (const auto placemark = geodata_cast<GeoDataPlacemark>(feature)) { 0030 // Select zoom level such that the placemark fits in a single tile 0031 int zoomLevel; 0032 qreal north, south, east, west; 0033 placemark->geometry()->latLonAltBox().boundaries(north, south, east, west); 0034 for (zoomLevel = maxZoomLevel; zoomLevel >= 0; --zoomLevel) { 0035 if (TileId::fromCoordinates(GeoDataCoordinates(west, north), zoomLevel) == 0036 TileId::fromCoordinates(GeoDataCoordinates(east, south), zoomLevel)) { 0037 break; 0038 } 0039 } 0040 TileId const key = TileId::fromCoordinates(GeoDataCoordinates(west, north), zoomLevel); 0041 m_items[key] << placemark; 0042 } else if (GeoDataRelation *relation = geodata_cast<GeoDataRelation>(feature)) { 0043 m_relations << relation; 0044 } else { 0045 Q_ASSERT(false && "only placemark variants are supported so far"); 0046 } 0047 } 0048 } 0049 0050 GeoDataDocument *VectorClipper::clipTo(const GeoDataLatLonBox &tileBoundary, int zoomLevel) 0051 { 0052 bool const filterSmallAreas = zoomLevel > 10 && zoomLevel < 17; 0053 GeoDataDocument* tile = new GeoDataDocument(); 0054 auto const clip = clipPath(tileBoundary, zoomLevel); 0055 GeoDataLinearRing ring; 0056 ring << GeoDataCoordinates(tileBoundary.west(), tileBoundary.north()); 0057 ring << GeoDataCoordinates(tileBoundary.east(), tileBoundary.north()); 0058 ring << GeoDataCoordinates(tileBoundary.east(), tileBoundary.south()); 0059 ring << GeoDataCoordinates(tileBoundary.west(), tileBoundary.south()); 0060 qreal const minArea = filterSmallAreas ? 0.01 * area(ring) : 0.0; 0061 QSet<qint64> osmIds; 0062 for (GeoDataPlacemark const * placemark: potentialIntersections(tileBoundary)) { 0063 GeoDataGeometry const * const geometry = placemark ? placemark->geometry() : nullptr; 0064 if (geometry && tileBoundary.intersects(geometry->latLonAltBox())) { 0065 if (geodata_cast<GeoDataPolygon>(geometry)) { 0066 clipPolygon(placemark, clip, minArea, tile, osmIds); 0067 } else if (geodata_cast<GeoDataLineString>(geometry)) { 0068 clipString<GeoDataLineString>(placemark, clip, minArea, tile, osmIds); 0069 } else if (geodata_cast<GeoDataLinearRing>(geometry)) { 0070 clipString<GeoDataLinearRing>(placemark, clip, minArea, tile, osmIds); 0071 } else if (const auto building = geodata_cast<GeoDataBuilding>(geometry)) { 0072 if (geodata_cast<GeoDataPolygon>(&static_cast<const GeoDataMultiGeometry*>(building->multiGeometry())->at(0))) { 0073 clipPolygon(placemark, clip, minArea, tile, osmIds); 0074 } else if (geodata_cast<GeoDataLinearRing>(&static_cast<const GeoDataMultiGeometry*>(building->multiGeometry())->at(0))) { 0075 clipString<GeoDataLinearRing>(placemark, clip, minArea, tile, osmIds); 0076 } 0077 } else { 0078 tile->append(placemark->clone()); 0079 osmIds << placemark->osmData().id(); 0080 } 0081 } 0082 } 0083 0084 for (auto relation: m_relations) { 0085 if (relation->containsAnyOf(osmIds)) { 0086 GeoDataRelation* multi = new GeoDataRelation; 0087 multi->osmData() = relation->osmData(); 0088 tile->append(multi); 0089 } 0090 } 0091 return tile; 0092 } 0093 0094 QVector<GeoDataPlacemark *> VectorClipper::potentialIntersections(const GeoDataLatLonBox &box) const 0095 { 0096 qreal north, south, east, west; 0097 box.boundaries(north, south, east, west); 0098 TileId const topLeft = TileId::fromCoordinates(GeoDataCoordinates(west, north), m_maxZoomLevel); 0099 TileId const bottomRight = TileId::fromCoordinates(GeoDataCoordinates(east, south), m_maxZoomLevel); 0100 QRect rect; 0101 rect.setCoords(topLeft.x(), topLeft.y(), bottomRight.x(), bottomRight.y()); 0102 0103 TileCoordsPyramid pyramid(0, m_maxZoomLevel); 0104 pyramid.setBottomLevelCoords(rect); 0105 QVector<GeoDataPlacemark *> result; 0106 for (int level = pyramid.topLevel(), maxLevel = pyramid.bottomLevel(); level <= maxLevel; ++level) { 0107 int x1, y1, x2, y2; 0108 pyramid.coords(level).getCoords(&x1, &y1, &x2, &y2); 0109 for (int x = x1; x <= x2; ++x) { 0110 for (int y = y1; y <= y2; ++y) { 0111 result << m_items.value(TileId(0, level, x, y)); 0112 } 0113 } 0114 } 0115 return result; 0116 } 0117 0118 GeoDataDocument *VectorClipper::clipTo(unsigned int zoomLevel, unsigned int tileX, unsigned int tileY) 0119 { 0120 const GeoDataLatLonBox tileBoundary = m_tileProjection.geoCoordinates(zoomLevel, tileX, tileY); 0121 0122 GeoDataDocument *tile = clipTo(tileBoundary, zoomLevel); 0123 QString tileName = QString("%1/%2/%3").arg(zoomLevel).arg(tileX).arg(tileY); 0124 tile->setName(tileName); 0125 0126 return tile; 0127 } 0128 0129 ClipperLib::Path VectorClipper::clipPath(const GeoDataLatLonBox &box, int zoomLevel) const 0130 { 0131 using namespace ClipperLib; 0132 Path path; 0133 int const steps = qMax(1, 22 - 2 * zoomLevel); 0134 double x = box.west() * s_pointScale; 0135 double const horizontalStep = (box.east() * s_pointScale - x) / steps; 0136 double y = box.north() * s_pointScale; 0137 double const verticalStep = (box.south() * s_pointScale - y) / steps; 0138 for (int i=0; i<steps; ++i) { 0139 path << IntPoint(qRound64(x), qRound64(y)); 0140 x += horizontalStep; 0141 } 0142 path << IntPoint(qRound64(box.east() * s_pointScale), qRound64(box.north() * s_pointScale)); 0143 for (int i=0; i<steps; ++i) { 0144 path << IntPoint(qRound64(x), qRound64(y)); 0145 y += verticalStep; 0146 } 0147 path << IntPoint(qRound64(box.east() * s_pointScale), qRound64(box.south() * s_pointScale)); 0148 for (int i=0; i<steps; ++i) { 0149 path << IntPoint(qRound64(x), qRound64(y)); 0150 x -= horizontalStep; 0151 } 0152 path << IntPoint(qRound64(box.west() * s_pointScale), qRound64(box.south() * s_pointScale)); 0153 for (int i=0; i<steps; ++i) { 0154 path << IntPoint(qRound64(x), qRound64(y)); 0155 y -= verticalStep; 0156 } 0157 path << IntPoint(qRound64(box.west() * s_pointScale), qRound64(box.north() * s_pointScale)); 0158 return path; 0159 } 0160 0161 bool VectorClipper::canBeArea(GeoDataPlacemark::GeoDataVisualCategory visualCategory) 0162 { 0163 if (visualCategory >= GeoDataPlacemark::HighwaySteps && visualCategory <= GeoDataPlacemark::HighwayMotorway) { 0164 return false; 0165 } 0166 if (visualCategory >= GeoDataPlacemark::RailwayRail && visualCategory <= GeoDataPlacemark::RailwayFunicular) { 0167 return false; 0168 } 0169 if (visualCategory >= GeoDataPlacemark::AdminLevel1 && visualCategory <= GeoDataPlacemark::AdminLevel11) { 0170 return false; 0171 } 0172 0173 if (visualCategory == GeoDataPlacemark::BoundaryMaritime || visualCategory == GeoDataPlacemark::InternationalDateLine) { 0174 return false; 0175 } 0176 0177 return true; 0178 } 0179 0180 qreal VectorClipper::area(const GeoDataLinearRing &ring) 0181 { 0182 int const n = ring.size(); 0183 qreal area = 0; 0184 if (n<3) { 0185 return area; 0186 } 0187 for (int i = 1; i < n; ++i ){ 0188 area += (ring[i].longitude() - ring[i-1].longitude() ) * ( ring[i].latitude() + ring[i-1].latitude()); 0189 } 0190 area += (ring[0].longitude() - ring[n-1].longitude() ) * (ring[0].latitude() + ring[n-1].latitude()); 0191 qreal const result = EARTH_RADIUS * EARTH_RADIUS * qAbs(area * 0.5); 0192 return result; 0193 } 0194 0195 void VectorClipper::getBounds(const ClipperLib::Path &path, ClipperLib::cInt &minX, ClipperLib::cInt &maxX, ClipperLib::cInt &minY, ClipperLib::cInt &maxY) const 0196 { 0197 Q_ASSERT(!path.empty()); 0198 minX = path[0].X; 0199 maxX = minX; 0200 minY = path[0].Y; 0201 maxY = minY; 0202 for (auto const & point: path) { 0203 if (point.X < minX) { 0204 minX = point.X; 0205 } else if (point.X > maxX) { 0206 maxX = point.X; 0207 } 0208 if (point.Y < minY) { 0209 minY = point.Y; 0210 } else if (point.Y > maxY) { 0211 maxY = point.Y; 0212 } 0213 } 0214 } 0215 0216 void VectorClipper::clipPolygon(const GeoDataPlacemark *placemark, const ClipperLib::Path &tileBoundary, qreal minArea, 0217 GeoDataDocument *document, QSet<qint64> &osmIds) 0218 { 0219 bool isBuilding = false; 0220 GeoDataPolygon* polygon; 0221 std::unique_ptr<GeoDataPlacemark> copyPlacemark; 0222 if (const auto building = geodata_cast<GeoDataBuilding>(placemark->geometry())) { 0223 polygon = geodata_cast<GeoDataPolygon>(&static_cast<GeoDataMultiGeometry*>(building->multiGeometry())->at(0)); 0224 isBuilding = true; 0225 } else { 0226 copyPlacemark.reset(new GeoDataPlacemark(*placemark)); 0227 polygon = geodata_cast<GeoDataPolygon>(copyPlacemark->geometry()); 0228 } 0229 0230 if (minArea > 0.0 && area(polygon->outerBoundary()) < minArea) { 0231 return; 0232 } 0233 using namespace ClipperLib; 0234 Path path; 0235 QHash<std::pair<cInt, cInt>, const GeoDataCoordinates*> coordMap; 0236 for(auto const & node: qAsConst(polygon)->outerBoundary()) { 0237 auto p = coordinateToPoint(node); 0238 coordMap.insert(std::make_pair(p.X, p.Y), &node); 0239 path.push_back(std::move(p)); 0240 } 0241 0242 cInt minX, maxX, minY, maxY; 0243 getBounds(tileBoundary, minX, maxX, minY, maxY); 0244 Clipper clipper; 0245 clipper.PreserveCollinear(true); 0246 clipper.AddPath(tileBoundary, ptClip, true); 0247 clipper.AddPath(path, ptSubject, true); 0248 Paths paths; 0249 clipper.Execute(ctIntersection, paths); 0250 for(const auto &path: paths) { 0251 GeoDataPlacemark* newPlacemark = new GeoDataPlacemark; 0252 newPlacemark->setVisible(placemark->isVisible()); 0253 newPlacemark->setVisualCategory(placemark->visualCategory()); 0254 GeoDataLinearRing outerRing; 0255 OsmPlacemarkData const & placemarkOsmData = placemark->osmData(); 0256 OsmPlacemarkData & newPlacemarkOsmData = newPlacemark->osmData(); 0257 int index = -1; 0258 OsmPlacemarkData const & outerRingOsmData = placemarkOsmData.memberReference(index); 0259 OsmPlacemarkData & newOuterRingOsmData = newPlacemarkOsmData.memberReference(index); 0260 pathToRing(path, &outerRing, outerRingOsmData, newOuterRingOsmData, coordMap); 0261 0262 GeoDataPolygon* newPolygon = new GeoDataPolygon; 0263 newPolygon->setOuterBoundary(outerRing); 0264 if (isBuilding) { 0265 const auto building = geodata_cast<GeoDataBuilding>(placemark->geometry()); 0266 GeoDataBuilding* newBuilding = new GeoDataBuilding(*building); 0267 newBuilding->multiGeometry()->clear(); 0268 newBuilding->multiGeometry()->append(newPolygon); 0269 newPlacemark->setGeometry(newBuilding); 0270 } else { 0271 newPlacemark->setGeometry(newPolygon); 0272 } 0273 if (placemarkOsmData.id() > 0) { 0274 newPlacemarkOsmData.addTag(QStringLiteral("mx:oid"), QString::number(placemarkOsmData.id())); 0275 } 0276 copyTags(placemarkOsmData, newPlacemarkOsmData); 0277 copyTags(outerRingOsmData, newOuterRingOsmData); 0278 if (outerRingOsmData.id() > 0) { 0279 newOuterRingOsmData.addTag(QStringLiteral("mx:oid"), QString::number(outerRingOsmData.id())); 0280 osmIds.insert(outerRingOsmData.id()); 0281 } 0282 0283 auto const & innerBoundaries = qAsConst(polygon)->innerBoundaries(); 0284 for (index = 0; index < innerBoundaries.size(); ++index) { 0285 auto const & innerBoundary = innerBoundaries.at(index); 0286 if (minArea > 0.0 && area(innerBoundary) < minArea) { 0287 continue; 0288 } 0289 0290 auto const & innerRingOsmData = placemarkOsmData.memberReference(index); 0291 clipper.Clear(); 0292 clipper.AddPath(path, ptClip, true); 0293 Path innerPath; 0294 coordMap.clear(); 0295 for(auto const & node: innerBoundary) { 0296 auto p = coordinateToPoint(node); 0297 coordMap.insert(std::make_pair(p.X, p.Y), &node); 0298 innerPath.push_back(std::move(p)); 0299 } 0300 clipper.AddPath(innerPath, ptSubject, true); 0301 Paths innerPaths; 0302 clipper.Execute(ctIntersection, innerPaths); 0303 for(auto const &innerPath: innerPaths) { 0304 int const newIndex = newPolygon->innerBoundaries().size(); 0305 auto & newInnerRingOsmData = newPlacemarkOsmData.memberReference(newIndex); 0306 GeoDataLinearRing innerRing; 0307 pathToRing(innerPath, &innerRing, innerRingOsmData, newInnerRingOsmData, coordMap); 0308 newPolygon->appendInnerBoundary(innerRing); 0309 if (innerRingOsmData.id() > 0) { 0310 newInnerRingOsmData.addTag(QStringLiteral("mx:oid"), QString::number(innerRingOsmData.id())); 0311 osmIds.insert(innerRingOsmData.id()); 0312 } 0313 copyTags(innerRingOsmData, newInnerRingOsmData); 0314 } 0315 } 0316 0317 OsmObjectManager::initializeOsmData(newPlacemark); 0318 document->append(newPlacemark); 0319 osmIds << placemark->osmData().id(); 0320 } 0321 } 0322 0323 void VectorClipper::copyTags(const GeoDataPlacemark &source, GeoDataPlacemark &target) const 0324 { 0325 copyTags(source.osmData(), target.osmData()); 0326 } 0327 0328 void VectorClipper::copyTags(const OsmPlacemarkData &originalPlacemarkData, OsmPlacemarkData &targetOsmData) const 0329 { 0330 for (auto iter=originalPlacemarkData.tagsBegin(), end=originalPlacemarkData.tagsEnd(); iter != end; ++iter) { 0331 targetOsmData.addTag(iter.key(), iter.value()); 0332 } 0333 } 0334 0335 }