44 #include <ogrsf_frmts.h>
46 #include <gdal_priv.h>
49 #ifdef CHECK_MEMORY_LEAKS
51 #endif // CHECK_MEMORY_LEAKS
64 myRTree(&
Triangle::addSelf), myRaster(0) {
88 WRITE_WARNING(
"Cannot supply height since no height data was loaded");
98 corners.push_back(
Position(floor(normX) + 0.5, floor(normY) + 0.5,
myRaster[(
int)normY * xSize + (
int)normX]));
99 if (normX - floor(normX) > 0.5) {
100 corners.push_back(
Position(floor(normX) + 1.5, floor(normY) + 0.5,
myRaster[(
int)normY * xSize + (
int)normX + 1]));
102 corners.push_back(
Position(floor(normX) - 0.5, floor(normY) + 0.5,
myRaster[(
int)normY * xSize + (
int)normX - 1]));
104 if (normY - floor(normY) > 0.5) {
105 corners.push_back(
Position(floor(normX) + 0.5, floor(normY) + 1.5,
myRaster[((
int)normY + 1) * xSize + (
int)normX]));
107 corners.push_back(
Position(floor(normX) + 0.5, floor(normY) - 0.5,
myRaster[((
int)normY - 1) * xSize + (
int)normX]));
111 if (result > -1e5 && result < 1e5) {
118 minB[0] = (float)geo.
x() - 0.00001f;
119 minB[1] = (float)geo.
y() - 0.00001f;
120 maxB[0] = (float)geo.
x() + 0.00001f;
121 maxB[1] = (float)geo.
y() + 0.00001f;
123 int hits =
myRTree.Search(minB, maxB, queryResult);
125 assert(hits == (
int)result.size());
128 for (Triangles::iterator it = result.begin(); it != result.end(); it++) {
131 return triangle->
getZ(geo);
144 const float cmin[2] = {(float) b.
xmin(), (float) b.
ymin()};
145 const float cmax[2] = {(float) b.
xmax(), (float) b.
ymax()};
146 myRTree.Insert(cmin, cmax, triangle);
152 if (oc.
isSet(
"heightmap.geotiff")) {
154 std::vector<std::string> files = oc.
getStringVector(
"heightmap.geotiff");
155 for (std::vector<std::string>::const_iterator file = files.begin(); file != files.end(); ++file) {
159 " done (parsed " +
toString(numFeatures) +
163 if (oc.
isSet(
"heightmap.shapefiles")) {
165 std::vector<std::string> files = oc.
getStringVector(
"heightmap.shapefiles");
166 for (std::vector<std::string>::const_iterator file = files.begin(); file != files.end(); ++file) {
170 " done (parsed " +
toString(numFeatures) +
180 #if GDAL_VERSION_MAJOR < 2
182 OGRDataSource* ds = OGRSFDriverRegistrar::Open(file.c_str(), FALSE);
185 GDALDataset* ds = (GDALDataset*)GDALOpenEx(file.c_str(), GDAL_OF_VECTOR | GA_ReadOnly, NULL, NULL, NULL);
188 throw ProcessError(
"Could not open shape file '" + file +
"'.");
192 OGRLayer* layer = ds->GetLayer(0);
193 layer->ResetReading();
197 OGRSpatialReference* sr_src = layer->GetSpatialRef();
198 OGRSpatialReference sr_dest;
199 sr_dest.SetWellKnownGeogCS(
"WGS84");
200 OGRCoordinateTransformation* toWGS84 = OGRCreateCoordinateTransformation(sr_src, &sr_dest);
202 WRITE_WARNING(
"Could not create geocoordinates converter; check whether proj.4 is installed.");
207 layer->ResetReading();
208 while ((feature = layer->GetNextFeature()) != NULL) {
209 OGRGeometry* geom = feature->GetGeometryRef();
213 assert(std::string(geom->getGeometryName()) == std::string(
"POLYGON"));
215 geom->transform(toWGS84);
216 OGRLinearRing* cgeom = ((OGRPolygon*) geom)->getExteriorRing();
218 assert(cgeom->getNumPoints() == 4);
220 for (
int j = 0; j < 3; j++) {
222 corners.push_back(pos);
259 OGRFeature::DestroyFeature(feature);
261 #if GDAL_VERSION_MAJOR < 2
262 OGRDataSource::DestroyDataSource(ds);
266 OCTDestroyCoordinateTransformation(toWGS84);
270 WRITE_ERROR(
"Cannot load shape file since SUMO was compiled without GDAL support.");
280 GDALDataset* poDataset = (GDALDataset*)GDALOpen(file.c_str(), GA_ReadOnly);
281 if (poDataset == 0) {
285 const int xSize = poDataset->GetRasterXSize();
286 const int ySize = poDataset->GetRasterYSize();
287 double adfGeoTransform[6];
288 if (poDataset->GetGeoTransform(adfGeoTransform) == CE_None) {
289 Position topLeft(adfGeoTransform[0], adfGeoTransform[3]);
294 myBoundary.
add(topLeft.
x() + horizontalSize, topLeft.
y() + verticalSize);
296 WRITE_ERROR(
"Could not parse geo information from " + file +
".");
299 const int picSize = xSize * ySize;
300 myRaster = (int16_t*)CPLMalloc(
sizeof(int16_t) * picSize);
301 for (
int i = 1; i <= poDataset->GetRasterCount(); i++) {
302 GDALRasterBand* poBand = poDataset->GetRasterBand(i);
303 if (poBand->GetColorInterpretation() != GCI_GrayIndex) {
304 WRITE_ERROR(
"Unknown color band in " + file +
".");
308 if (poBand->GetRasterDataType() != GDT_Int16) {
313 assert(xSize == poBand->GetXSize() && ySize == poBand->GetYSize());
314 if (poBand->RasterIO(GF_Read, 0, 0, xSize, ySize,
myRaster, xSize, ySize, GDT_Int16, 0, 0) == CE_Failure) {
320 GDALClose(poDataset);
323 WRITE_ERROR(
"Cannot load GeoTIFF file since SUMO was compiled without GDAL support.");
363 return myCorners.around(pos);
381 Position side1 = myCorners[1] - myCorners[0];
382 Position side2 = myCorners[2] - myCorners[0];
void sub(SUMOReal dx, SUMOReal dy)
Substracts the given position from this one.
std::vector< std::string > getStringVector(const std::string &name) const
Returns the list of string-vector-value of the named option (only for Option_String) ...
SUMOReal getZ(const Position &geo) const
returns the projection of the give geoCoordinate (WGS84) onto triangle plane
bool ready() const
returns whether the NBHeightMapper has data
bool around(const Position &p, SUMOReal offset=0) const
Returns whether the boundary contains the given coordinate.
Position normalVector() const
returns the normal vector for this triangles plane
int loadShapeFile(const std::string &file)
load height data from Arcgis-shape file and returns the number of parsed features ...
SUMOReal getZ(const Position &geo) const
returns height for the given geo coordinate (WGS84)
SUMOReal ymin() const
Returns minimum y-coordinate.
SUMOReal xmin() const
Returns minimum x-coordinate.
Triangle(const PositionVector &corners)
void addTriangle(PositionVector corners)
adds one triangles worth of height data
void clearData()
clears loaded data
static NBHeightMapper Singleton
the singleton instance
NBHeightMapper()
private constructor and destructor (Singleton)
SUMOReal x() const
Returns the x-position.
#define UNUSED_PARAMETER(x)
SUMOReal xmax() const
Returns maximum x-coordinate.
A class that stores a 2D geometrical boundary.
#define WRITE_WARNING(msg)
int loadTiff(const std::string &file)
load height data from GeoTIFF file and returns the number of non void pixels
static const NBHeightMapper & get()
return the singleton instance (maybe 0)
class for cirumventing the const-restriction of RTree::Search-context
static void loadIfSet(OptionsCont &oc)
loads heigh map data if any loading options are set
TRIANGLE_RTREE_QUAL myRTree
The RTree for spatial queries.
Position crossProduct(const Position &pos)
returns the cross product between this point and the second one
A point in 2D or 3D with translation and scaling methods.
std::vector< const Triangle * > Triangles
#define PROGRESS_BEGIN_MESSAGE(msg)
static MsgHandler * getMessageInstance()
Returns the instance to add normal messages to.
int16_t * myRaster
raster height information in m
bool contains(const Position &pos) const
checks whether pos lies within triangle (only checks x,y)
std::string toString(const T &t, std::streamsize accuracy=OUTPUT_ACCURACY)
Position mySizeOfPixel
dimensions of one pixel in raster data
Boundary myBoundary
convex boundary of all known triangles;
void reset()
Resets the boundary.
void add(SUMOReal x, SUMOReal y, SUMOReal z=0)
Makes the boundary include the given coordinate.
SUMOReal dotProduct(const Position &pos)
returns the dot product (scalar product) between this point and the second one
SUMOReal y() const
Returns the y-position.
PositionVector myCorners
the corners of the triangle
A storage for options typed value containers)
void set(SUMOReal x, SUMOReal y)
void addSelf(const QueryResult &queryResult) const
callback for RTree search
SUMOReal ymax() const
Returns maximum y-coordinate.
const Boundary & getBoundary()
returns the convex boundary of all known triangles
Boundary getBoxBoundary() const
Returns a boundary enclosing this list of lines.
Set z-values for all network positions based on data from a height map.
bool isSet(const std::string &name, bool failOnNonExistant=true) const
Returns the information whether the named option is set.
void endProcessMsg(std::string msg)
Ends a process information.