From 2b22e7b361c12624dd0f1f25a88f6399107469d1 Mon Sep 17 00:00:00 2001 From: Dakota Benjamin Date: Fri, 26 Jan 2018 14:38:26 -0500 Subject: [PATCH] Refactor georeferencing to implement sfm transformation --- SuperBuild/cmake/External-OpenSfM.cmake | 4 +- modules/odm_georef/src/Georef.cpp | 336 +++++++++---------- modules/odm_georef/src/Georef.hpp | 5 + modules/odm_orthophoto/src/OdmOrthoPhoto.cpp | 90 ++++- modules/odm_orthophoto/src/OdmOrthoPhoto.hpp | 9 + opendm/config.py | 4 + opendm/io.py | 6 + opendm/types.py | 122 ++++--- scripts/dataset.py | 50 ++- scripts/mvstex.py | 12 +- scripts/odm_app.py | 9 +- scripts/odm_georeferencing.py | 89 ++--- scripts/odm_meshing.py | 7 +- scripts/odm_orthophoto.py | 14 +- scripts/{opensfm.py => run_opensfm.py} | 20 +- 15 files changed, 454 insertions(+), 323 deletions(-) rename scripts/{opensfm.py => run_opensfm.py} (92%) diff --git a/SuperBuild/cmake/External-OpenSfM.cmake b/SuperBuild/cmake/External-OpenSfM.cmake index 9fa30f42..874383c4 100644 --- a/SuperBuild/cmake/External-OpenSfM.cmake +++ b/SuperBuild/cmake/External-OpenSfM.cmake @@ -8,8 +8,8 @@ ExternalProject_Add(${_proj_name} STAMP_DIR ${_SB_BINARY_DIR}/stamp #--Download step-------------- DOWNLOAD_DIR ${SB_DOWNLOAD_DIR} - URL https://github.com/mapillary/OpenSfM/archive/93be3a1bfe46482345ddd57bc1f3a62f63169b86.zip - URL_MD5 2b310420a5c7c2297294a39183fb8b1a + URL https://github.com/mapillary/OpenSfM/archive/bcffb387b1d41a0c82919c92fd38b07d1c2bed90.zip + # URL_MD5 72587e6b7011e746aa4208d785e4b49a #--Update/Patch step---------- UPDATE_COMMAND "" #--Configure step------------- diff --git a/modules/odm_georef/src/Georef.cpp b/modules/odm_georef/src/Georef.cpp index 676c0587..b77cbeb1 100644 --- a/modules/odm_georef/src/Georef.cpp +++ b/modules/odm_georef/src/Georef.cpp @@ -236,6 +236,7 @@ Georef::Georef() : log_(false) outputCoordFilename_ = ""; inputObjFilename_ = ""; outputObjFilename_ = ""; + transformFilename_ = ""; exportCoordinateFile_ = false; exportGeorefSystem_ = false; } @@ -416,6 +417,17 @@ void Georef::parseArguments(int argc, char *argv[]) log_ << "Reading GCPs from: " << gcpFilename_ << "\n"; gcpFileSpecified = true; } + else if(argument == "-inputTransformFile" && argIndex < argc) + { + argIndex++; + if (argIndex >= argc) + { + throw GeorefException("Argument '" + argument + "' expects 1 more input following it, but no more inputs were provided."); + } + transformFilename_ = std::string(argv[argIndex]); + log_ << "Reading transform file from: " << gcpFilename_ << "\n"; + useTransform_ = true; + } else if(argument == "-imagesListPath" && argIndex < argc) { argIndex++; @@ -558,6 +570,9 @@ void Georef::printHelp() log_ << "\"-inputPointCloudFile \" (optional)" << "\n"; log_ << "\"Input ply file that must contain a point cloud.\n\n"; + log_ << "\"-inputTransformFile \" (optional)" << "\n"; + log_ << "\"Input transform file that is a 4x4 matrix of transform values.\n\n"; + log_ << "\"-imagesListPath \" (mandatory if using ground control points)\n"; log_ << "Path to the list containing the image names used in the bundle.out file.\n\n"; @@ -624,6 +639,10 @@ void Georef::createGeoreferencedModel() { createGeoreferencedModelFromGCPData(); } + else if (useTransform_) + { + createGeoreferencedModelFromSFM(); + } else { createGeoreferencedModelFromExifData(); @@ -911,107 +930,7 @@ void Georef::performGeoreferencingWithGCP() printFinalTransform(transFinal.transform_); // The transform used to transform model into the georeferenced system. - Eigen::Transform transform; - - transform(0, 0) = static_cast(transFinal.transform_.r1c1_); - transform(1, 0) = static_cast(transFinal.transform_.r2c1_); - transform(2, 0) = static_cast(transFinal.transform_.r3c1_); - transform(3, 0) = static_cast(transFinal.transform_.r4c1_); - - transform(0, 1) = static_cast(transFinal.transform_.r1c2_); - transform(1, 1) = static_cast(transFinal.transform_.r2c2_); - transform(2, 1) = static_cast(transFinal.transform_.r3c2_); - transform(3, 1) = static_cast(transFinal.transform_.r4c2_); - - transform(0, 2) = static_cast(transFinal.transform_.r1c3_); - transform(1, 2) = static_cast(transFinal.transform_.r2c3_); - transform(2, 2) = static_cast(transFinal.transform_.r3c3_); - transform(3, 2) = static_cast(transFinal.transform_.r4c3_); - - transform(0, 3) = static_cast(transFinal.transform_.r1c4_); - transform(1, 3) = static_cast(transFinal.transform_.r2c4_); - transform(2, 3) = static_cast(transFinal.transform_.r3c4_); - transform(3, 3) = static_cast(transFinal.transform_.r4c4_); - - log_ << '\n'; - log_ << "Applying transform to mesh...\n"; - // Move the mesh into position. - pcl::transformPointCloud(*meshCloud, *meshCloud, transform); - log_ << ".. mesh transformed.\n"; - - // Update the mesh. - pcl::toPCLPointCloud2 (*meshCloud, mesh.cloud); - - // Iterate over each part of the mesh (one per material), to make texture file paths relative the .mtl file. - for(size_t t = 0; t < mesh.tex_materials.size(); ++t) - { - // The material of the current submesh. - pcl::TexMaterial& material = mesh.tex_materials[t]; - - size_t find = material.tex_file.find_last_of("/\\"); - if(std::string::npos != find) - { - material.tex_file = material.tex_file.substr(find + 1); - } - } - - log_ << '\n'; - if (saveOBJFile(outputObjFilename_, mesh, 8) == -1) - { - throw GeorefException("Error when saving model:\n" + outputObjFilename_ + "\n"); - } - else - { - log_ << "Successfully saved model.\n"; - } - - if(georeferencePointCloud_) - { - //pcl::PointCloud2::Ptr pointCloud; - pcl::PointCloud::Ptr pointCloud(new pcl::PointCloud()); - if(pcl::io::loadPLYFile (inputPointCloudFilename_.c_str(), *pointCloud.get()) == -1) { - throw GeorefException("Error when reading point cloud:\n" + inputPointCloudFilename_ + "\n"); - } - else - { - log_ << "Successfully loaded " << pointCloud->size() << " points with corresponding normals from file.\n"; - } - log_ << '\n'; - log_ << "Applying transform to point cloud...\n"; - pcl::transformPointCloud(*pointCloud, *pointCloud, transform); - log_ << ".. point cloud transformed.\n"; - - pcl::PLYWriter plyWriter; - - log_ << '\n'; - log_ << "Saving point cloud file to \'" << outputPointCloudFilename_ << "\'...\n"; - //pcl::io::savePLYFileASCII(outputPointCloudFilename_.c_str(), *pointCloud.get()); - plyWriter.write(outputPointCloudFilename_.c_str(), *pointCloud.get(), false, false); - log_ << ".. point cloud file saved.\n"; - } - - if(exportCoordinateFile_) - { - log_ << '\n'; - log_ << "Saving georeferenced camera positions to "; - log_ << outputCoordFilename_; - log_<< "\n"; - std::ofstream coordStream(outputCoordFilename_.c_str()); - coordStream << georefSystem_.system_ <(georefSystem_.eastingOffset_) << " " << static_cast(georefSystem_.northingOffset_) << std::endl; - for(size_t cameraIndex = 0; cameraIndex < cameras_.size(); ++cameraIndex) - { - Vec3 globalCameraPosition = (transFinal.transform_)*(cameras_[cameraIndex].getPos()); - coordStream << globalCameraPosition.x_ << " " << globalCameraPosition.y_ << " " << globalCameraPosition.z_ << std::endl; - } - coordStream.close(); - log_ << "...coordinate file saved.\n"; - } - - if(exportGeorefSystem_) - { - printGeorefSystem(); - } + performFinalTransform(transFinal.transform_, mesh, meshCloud); } void Georef::createGeoreferencedModelFromGCPData() @@ -1097,88 +1016,67 @@ void Georef::createGeoreferencedModelFromExifData() log_ << transFinal.transform_ << '\n'; printFinalTransform(transFinal.transform_); - - // The transform used to move the chosen area into the ortho photo. - Eigen::Transform transform; - - transform(0, 0) = static_cast(transFinal.transform_.r1c1_); transform(1, 0) = static_cast(transFinal.transform_.r2c1_); - transform(2, 0) = static_cast(transFinal.transform_.r3c1_); transform(3, 0) = static_cast(transFinal.transform_.r4c1_); - - transform(0, 1) = static_cast(transFinal.transform_.r1c2_); transform(1, 1) = static_cast(transFinal.transform_.r2c2_); - transform(2, 1) = static_cast(transFinal.transform_.r3c2_); transform(3, 1) = static_cast(transFinal.transform_.r4c2_); - - transform(0, 2) = static_cast(transFinal.transform_.r1c3_); transform(1, 2) = static_cast(transFinal.transform_.r2c3_); - transform(2, 2) = static_cast(transFinal.transform_.r3c3_); transform(3, 2) = static_cast(transFinal.transform_.r4c3_); - - transform(0, 3) = static_cast(transFinal.transform_.r1c4_); transform(1, 3) = static_cast(transFinal.transform_.r2c4_); - transform(2, 3) = static_cast(transFinal.transform_.r3c4_); transform(3, 3) = static_cast(transFinal.transform_.r4c4_); - + log_ << '\n'; log_ << "Reading mesh file...\n"; pcl::TextureMesh mesh; loadObjFile(inputObjFilename_, mesh); log_ << ".. mesh file read.\n"; - + // Contains the vertices of the mesh. pcl::PointCloud::Ptr meshCloud (new pcl::PointCloud); pcl::fromPCLPointCloud2 (mesh.cloud, *meshCloud); - + + performFinalTransform(transFinal.transform_, mesh, meshCloud); +} + +void Georef::createGeoreferencedModelFromSFM() +{ + + Mat4 transform; + // get transform in correct format + // Read elements from transform file generated by opensfm + std::ifstream transStream(transformFilename_.c_str()); + if (!transStream.good()) + { + throw GeorefException("Failed opening coordinate file " + transformFilename_ + " for reading. " + '\n'); + } + + std::string transString; + { + std::getline(transStream, transString); + std::stringstream l1(transString); + l1 >> transform.r1c1_ >> transform.r1c2_ >> transform.r1c3_ >> transform.r1c4_; + + std::getline(transStream, transString); + std::stringstream l2(transString); + l2 >> transform.r2c1_ >> transform.r2c2_ >> transform.r2c3_ >> transform.r2c4_; + + std::getline(transStream, transString); + std::stringstream l3(transString); + l3 >> transform.r3c1_ >> transform.r3c2_ >> transform.r3c3_ >> transform.r3c4_; + + std::getline(transStream, transString); + std::stringstream l4(transString); + l4 >> transform.r4c1_ >> transform.r4c2_ >> transform.r4c3_ >> transform.r4c4_; + } + + // load mesh + printFinalTransform(transform); + log_ << '\n'; - log_ << "Applying transform to mesh...\n"; - // Move the mesh into position. - pcl::transformPointCloud(*meshCloud, *meshCloud, transform); - log_ << ".. mesh transformed.\n"; - - // Update the mesh. - pcl::toPCLPointCloud2 (*meshCloud, mesh.cloud); - - // Iterate over each part of the mesh (one per material), to make texture file paths relative the .mtl file. - for(size_t t = 0; t < mesh.tex_materials.size(); ++t) - { - // The material of the current submesh. - pcl::TexMaterial& material = mesh.tex_materials[t]; - - size_t find = material.tex_file.find_last_of("/\\"); - if(std::string::npos != find) - { - material.tex_file = material.tex_file.substr(find + 1); - } - } - - log_ << '\n'; - log_ << "Saving mesh file to \'" << outputObjFilename_ << "\'...\n"; - saveOBJFile(outputObjFilename_, mesh, 8); - log_ << ".. mesh file saved.\n"; - - if(georeferencePointCloud_) - { - //pcl::PointCloud2::Ptr pointCloud; - pcl::PointCloud::Ptr pointCloud(new pcl::PointCloud()); - if(pcl::io::loadPLYFile (inputPointCloudFilename_.c_str(), *pointCloud.get()) == -1) { - throw GeorefException("Error when reading point cloud:\n" + inputPointCloudFilename_ + "\n"); - } - else - { - log_ << "Successfully loaded " << pointCloud->size() << " points with corresponding normals from file.\n"; - } - log_ << '\n'; - log_ << "Applying transform to point cloud...\n"; - pcl::transformPointCloud(*pointCloud, *pointCloud, transform); - log_ << ".. point cloud transformed.\n"; + log_ << "Reading mesh file...\n"; + pcl::TextureMesh mesh; + loadObjFile(inputObjFilename_, mesh); + log_ << ".. mesh file read.\n"; - pcl::PLYWriter plyWriter; + // Contains the vertices of the mesh. + pcl::PointCloud::Ptr meshCloud (new pcl::PointCloud); + pcl::fromPCLPointCloud2 (mesh.cloud, *meshCloud); - log_ << '\n'; - log_ << "Saving point cloud file to \'" << outputPointCloudFilename_ << "\'...\n"; - //pcl::io::savePLYFileASCII(outputPointCloudFilename_.c_str(), *pointCloud.get()); - plyWriter.write(outputPointCloudFilename_.c_str(), *pointCloud.get(), false, false); - log_ << ".. point cloud file saved.\n"; - } + performFinalTransform(transform, mesh, meshCloud); - if(exportGeorefSystem_) - { - printGeorefSystem(); - } + // performFinalTransform(transFinal, mesh, meshCloud); } void Georef::chooseBestGCPTriplet(size_t &gcp0, size_t &gcp1, size_t &gcp2) @@ -1377,6 +1275,102 @@ void Georef::printFinalTransform(Mat4 transform) } +void Georef::performFinalTransform(Mat4 &transMat, pcl::TextureMesh &mesh, pcl::PointCloud::Ptr &meshCloud) +{ + Eigen::Transform transform; + + transform(0, 0) = static_cast(transMat.r1c1_); + transform(1, 0) = static_cast(transMat.r2c1_); + transform(2, 0) = static_cast(transMat.r3c1_); + transform(3, 0) = static_cast(transMat.r4c1_); + + transform(0, 1) = static_cast(transMat.r1c2_); + transform(1, 1) = static_cast(transMat.r2c2_); + transform(2, 1) = static_cast(transMat.r3c2_); + transform(3, 1) = static_cast(transMat.r4c2_); + + transform(0, 2) = static_cast(transMat.r1c3_); + transform(1, 2) = static_cast(transMat.r2c3_); + transform(2, 2) = static_cast(transMat.r3c3_); + transform(3, 2) = static_cast(transMat.r4c3_); + + transform(0, 3) = static_cast(transMat.r1c4_); + transform(1, 3) = static_cast(transMat.r2c4_); + transform(2, 3) = static_cast(transMat.r3c4_); + transform(3, 3) = static_cast(transMat.r4c4_); + + // Iterate over each part of the mesh (one per material), to make texture file paths relative the .mtl file. + for(size_t t = 0; t < mesh.tex_materials.size(); ++t) + { + // The material of the current submesh. + pcl::TexMaterial& material = mesh.tex_materials[t]; + + size_t find = material.tex_file.find_last_of("/\\"); + if(std::string::npos != find) + { + material.tex_file = material.tex_file.substr(find + 1); + } + } + + log_ << '\n'; + if (saveOBJFile(outputObjFilename_, mesh, 8) == -1) + { + throw GeorefException("Error when saving model:\n" + outputObjFilename_ + "\n"); + } + else + { + log_ << "Successfully saved model.\n"; + } + + if(georeferencePointCloud_) + { + //pcl::PointCloud2::Ptr pointCloud; + pcl::PointCloud::Ptr pointCloud(new pcl::PointCloud()); + if(pcl::io::loadPLYFile (inputPointCloudFilename_.c_str(), *pointCloud.get()) == -1) { + throw GeorefException("Error when reading point cloud:\n" + inputPointCloudFilename_ + "\n"); + } + else + { + log_ << "Successfully loaded " << pointCloud->size() << " points with corresponding normals from file.\n"; + } + log_ << '\n'; + log_ << "Applying transform to point cloud...\n"; + pcl::transformPointCloud(*pointCloud, *pointCloud, transform); + log_ << ".. point cloud transformed.\n"; + + pcl::PLYWriter plyWriter; + + log_ << '\n'; + log_ << "Saving point cloud file to \'" << outputPointCloudFilename_ << "\'...\n"; + //pcl::io::savePLYFileASCII(outputPointCloudFilename_.c_str(), *pointCloud.get()); + plyWriter.write(outputPointCloudFilename_.c_str(), *pointCloud.get(), false, false); + log_ << ".. point cloud file saved.\n"; + } + + if(exportCoordinateFile_) + { + log_ << '\n'; + log_ << "Saving georeferenced camera positions to "; + log_ << outputCoordFilename_; + log_<< "\n"; + std::ofstream coordStream(outputCoordFilename_.c_str()); + coordStream << georefSystem_.system_ <(georefSystem_.eastingOffset_) << " " << static_cast(georefSystem_.northingOffset_) << std::endl; + for(size_t cameraIndex = 0; cameraIndex < cameras_.size(); ++cameraIndex) + { + Vec3 globalCameraPosition = (transMat)*(cameras_[cameraIndex].getPos()); + coordStream << globalCameraPosition.x_ << " " << globalCameraPosition.y_ << " " << globalCameraPosition.z_ << std::endl; + } + coordStream.close(); + log_ << "...coordinate file saved.\n"; + } + + if(exportGeorefSystem_) + { + printGeorefSystem(); + } +} + bool Georef::loadObjFile(std::string inputFile, pcl::TextureMesh &mesh) { int data_type; diff --git a/modules/odm_georef/src/Georef.hpp b/modules/odm_georef/src/Georef.hpp index 17b09a6e..b299cbd5 100644 --- a/modules/odm_georef/src/Georef.hpp +++ b/modules/odm_georef/src/Georef.hpp @@ -264,6 +264,7 @@ private: std::string inputCoordFilename_; /**< The path to the cameras exif gps positions file. **/ std::string outputCoordFilename_; /**< The path to the cameras georeferenced gps positions file. **/ std::string gcpFilename_; /**< The path to the GCP file **/ + std::string transformFilename_; /**< The path to the input transform file **/ std::string imagesListPath_; /**< Path to the image list. **/ std::string imagesLocation_; /**< The folder containing the images in the image list. **/ std::string inputObjFilename_; /**< The path to the input mesh obj file. **/ @@ -276,6 +277,7 @@ private: bool exportCoordinateFile_; bool exportGeorefSystem_; bool useGCP_; /**< Check if GCP-file is present and use this to georeference the model. **/ + bool useTransform_; // double bundleResizedTo_; /**< The size used in the previous steps to calculate the camera focal_length. */ std::vector cameras_; /**< A vector of all cameras. **/ @@ -287,6 +289,9 @@ private: bool multiMaterial_; /**< True if the mesh has multiple materials. **/ std::vector companions_; /**< Materials (used by loadOBJFile). **/ + void performFinalTransform(Mat4 &transMat, pcl::TextureMesh &mesh, pcl::PointCloud::Ptr &meshCloud); + + void createGeoreferencedModelFromSFM(); }; /*! diff --git a/modules/odm_orthophoto/src/OdmOrthoPhoto.cpp b/modules/odm_orthophoto/src/OdmOrthoPhoto.cpp index 906f5464..5f3558f2 100644 --- a/modules/odm_orthophoto/src/OdmOrthoPhoto.cpp +++ b/modules/odm_orthophoto/src/OdmOrthoPhoto.cpp @@ -44,10 +44,13 @@ OdmOrthoPhoto::OdmOrthoPhoto() { inputFile_ = ""; inputGeoRefFile_ = ""; + inputTransformFile_ = ""; outputFile_ = "ortho.jpg"; logFile_ = "log.txt"; outputCornerFile_ = ""; + transformOverride_ = false; + resolution_ = 0.0f; boundaryDefined_ = false; @@ -242,6 +245,17 @@ void OdmOrthoPhoto::parseArguments(int argc, char *argv[]) outputCornerFile_ = std::string(argv[argIndex]); log_ << "Writing corners to: " << outputCornerFile_ << "\n"; } + else if(argument == "-inputTransformFile") + { + argIndex++; + if (argIndex >= argc) + { + throw OdmOrthoPhotoException("Argument '" + argument + "' expects 1 more input following it, but no more inputs were provided."); + } + inputTransformFile_ = std::string(argv[argIndex]); + transformOverride_ = true; + log_ << "Reading transformation matrix from: " << outputCornerFile_ << "\n"; + } else { printHelp(); @@ -655,28 +669,70 @@ void OdmOrthoPhoto::adjustBoundsForEntireModel(const pcl::TextureMesh &mesh) Eigen::Transform OdmOrthoPhoto::getROITransform(float xMin, float yMin) const { - // The transform used to move the chosen area into the ortho photo. + //Use transformation matrix if provided: + if(transformOverride_){ + return readTransform(inputTransformFile_); + } + else { + // The transform used to move the chosen area into the ortho photo. + Eigen::Transform transform; + + transform(0, 0) = resolution_; // x Scaling. + transform(1, 0) = 0.0f; + transform(2, 0) = 0.0f; + transform(3, 0) = 0.0f; + + transform(0, 1) = 0.0f; + transform(1, 1) = -resolution_; // y Scaling, mirrored for easier rendering. + transform(2, 1) = 0.0f; + transform(3, 1) = 0.0f; + + transform(0, 2) = 0.0f; + transform(1, 2) = 0.0f; + transform(2, 2) = 1.0f; + transform(3, 2) = 0.0f; + + transform(0, 3) = -xMin * resolution_; // x Translation + transform(1, 3) = -yMin * resolution_; // y Translation + transform(2, 3) = 0.0f; + transform(3, 3) = 1.0f; + + return transform; + } +} + +Eigen::Transform OdmOrthoPhoto::readTransform(std::string transformFile_) const +{ Eigen::Transform transform; - transform(0, 0) = resolution_; // x Scaling. - transform(1, 0) = 0.0f; - transform(2, 0) = 0.0f; - transform(3, 0) = 0.0f; + std::ifstream transStream(transformFile_.c_str()); + if (!transStream.good()) + { + throw OdmOrthoPhotoException("Failed opening coordinate file " + transformFile_ + " for reading. " + '\n'); + } - transform(0, 1) = 0.0f; - transform(1, 1) = -resolution_; // y Scaling, mirrored for easier rendering. - transform(2, 1) = 0.0f; - transform(3, 1) = 0.0f; + std::string transString; + { + std::getline(transStream, transString); + std::stringstream l1(transString); + l1 >> transform(0,0) >> transform(0,1) >> transform(0,2) >> transform(0,3); - transform(0, 2) = 0.0f; - transform(1, 2) = 0.0f; - transform(2, 2) = 1.0f; - transform(3, 2) = 0.0f; + std::getline(transStream, transString); + std::stringstream l2(transString); + l2 >> transform(1,0) >> transform(1,1) >> transform(1,2) >> transform(1,3); - transform(0, 3) = -xMin*resolution_; // x Translation - transform(1, 3) = -yMin*resolution_; // y Translation - transform(2, 3) = 0.0f; - transform(3, 3) = 1.0f; + std::getline(transStream, transString); + std::stringstream l3(transString); + l3 >> transform(2,0) >> transform(2,1) >> transform(2,2) >> transform(2,3); + + std::getline(transStream, transString); + std::stringstream l4(transString); + l4 >> transform(3,0) >> transform(3,1) >> transform(3,2) >> transform(3,3); + } + + // Don't do any rotation/shear + transform(0,1) = 0.0f; + transform(1,0) = 0.0f; return transform; } diff --git a/modules/odm_orthophoto/src/OdmOrthoPhoto.hpp b/modules/odm_orthophoto/src/OdmOrthoPhoto.hpp index 1b83f0c0..81a58201 100644 --- a/modules/odm_orthophoto/src/OdmOrthoPhoto.hpp +++ b/modules/odm_orthophoto/src/OdmOrthoPhoto.hpp @@ -116,6 +116,13 @@ private: * \brief Creates a transformation which aligns the area for the orthophoto. */ Eigen::Transform getROITransform(float xMin, float yMin) const; + + /*! + * \brief Reads a transformation matrix from a file + * @param transformFile_ + * @return + */ + Eigen::Transform readTransform(std::string transformFile_) const; /*! * \brief Renders a triangle into the ortho photo. @@ -193,12 +200,14 @@ private: std::string inputFile_; /**< Path to the textured mesh as an obj-file. */ std::string inputGeoRefFile_; /**< Path to the georeference system file. */ + std::string inputTransformFile_; std::string outputFile_; /**< Path to the destination file. */ std::string outputCornerFile_; /**< Path to the output corner file. */ std::string logFile_; /**< Path to the log file. */ float resolution_; /**< The number of pixels per meter in the ortho photo. */ + bool transformOverride_; bool boundaryDefined_; /**< True if the user has defined a boundary. */ WorldPoint worldPoint1_; /**< The first boundary point for the ortho photo, in world coordinates. */ diff --git a/opendm/config.py b/opendm/config.py index 5e7e91a3..b5ea9a4e 100644 --- a/opendm/config.py +++ b/opendm/config.py @@ -97,6 +97,10 @@ def config(): help=('Override the focal length information for the ' 'images')) + parser.add_argument('--proj', + metavar='', + help='Projection used to transform the model into geographic coordinates') + parser.add_argument('--force-ccd', metavar='', type=float, diff --git a/opendm/io.py b/opendm/io.py index fa566c44..6339b820 100644 --- a/opendm/io.py +++ b/opendm/io.py @@ -40,3 +40,9 @@ def copy(src, dst): if e.errno == errno.ENOTDIR: shutil.copy(src, dst) else: raise + + +# find a file in the root directory +def find(filename, folder): + for root, dirs, files in os.walk(folder): + return '/'.join((root, filename)) if filename in files else None diff --git a/opendm/types.py b/opendm/types.py index a981c93d..ed444d75 100644 --- a/opendm/types.py +++ b/opendm/types.py @@ -3,6 +3,7 @@ import pyexiv2 import re from fractions import Fraction from opensfm.exif import sensor_string +from pyproj import Proj import log import io @@ -128,9 +129,48 @@ class ODM_Photo: class ODM_Reconstruction(object): """docstring for ODMReconstruction""" - def __init__(self, arg): - super(ODMReconstruction, self).__init__() - self.arg = arg + def __init__(self, photos, projstring = None, coords_file = None): + self.photos = photos # list of ODM_Photos + self.projection = None # Projection system the whole project will be in + if projstring: + self.projection = self.set_projection(projstring) + self.georef = ODM_GeoRef(self.projection) + else: + self.projection = self.parse_coordinate_system(coords_file) + self.georef = ODM_GeoRef(self.projection) + + def parse_coordinate_system(self, _file): + """Write attributes to jobOptions from coord file""" + # check for coordinate file existence + if not io.file_exists(_file): + log.ODM_ERROR('Could not find file %s' % _file) + return + + with open(_file) as f: + # extract reference system and utm zone from first line. + # We will assume the following format: + # 'WGS84 UTM 17N' or 'WGS84 UTM 17N \n' + line = f.readline().rstrip() + log.ODM_DEBUG('Line: %s' % line) + ref = line.split(' ') + # match_wgs_utm = re.search('WGS84 UTM (\d{1,2})(N|S)', line, re.I) + if ref[0] == 'WGS84' and ref[1] == 'UTM': # match_wgs_utm: + datum = ref[0] + utm_pole = ref[2][len(ref[2]) - 1] + utm_zone = int(ref[2][:len(ref[2]) - 1]) + + return Proj(proj="utm", zone=utm_zone, datum=datum, no_defs=True) + elif '+proj' in line: + return Proj(line.strip('\'')) + else: + log.ODM_ERROR('Could not parse coordinates. Bad CRS supplied: %s' % line) + return + + def set_projection(self, projstring): + try: + return Proj(projstring) + except RuntimeError: + log.ODM_EXCEPTION('Could not set projection. Please use a proj4 string') class ODM_GCPoint(object): @@ -145,13 +185,15 @@ class ODM_GCPoint(object): class ODM_GeoRef(object): """docstring for ODMUtmZone""" - def __init__(self): + def __init__(self, projection): + self.projection = projection self.datum = 'WGS84' self.epsg = None self.utm_zone = 0 self.utm_pole = 'N' self.utm_east_offset = 0 self.utm_north_offset = 0 + self.transform = [] self.gcps = [] def calculate_EPSG(self, _utm_zone, _pole): @@ -164,6 +206,9 @@ class ODM_GeoRef(object): log.ODM_ERROR('Unknown pole format %s' % _pole) return + def calculate_EPSG(self, proj): + return proj + def coord_to_fractions(self, coord, refs): deg_dec = abs(float(coord)) deg = int(deg_dec) @@ -184,8 +229,8 @@ class ODM_GeoRef(object): def convert_to_las(self, _file, _file_out, json_file): - if not self.epsg: - log.ODM_ERROR('Empty EPSG: Could not convert to LAS') + if not self.projection.srs: + log.ODM_ERROR('Empty CRS: Could not convert to LAS') return kwargs = {'bin': context.pdal_path, @@ -193,7 +238,7 @@ class ODM_GeoRef(object): 'f_out': _file_out, 'east': self.utm_east_offset, 'north': self.utm_north_offset, - 'epsg': self.epsg, + 'srs': self.projection.srs, 'json': json_file} # create pipeline file transform.xml to enable transformation @@ -205,7 +250,7 @@ class ODM_GeoRef(object): ' "matrix":"1 0 0 {east} 0 1 0 {north} 0 0 1 0 0 0 0 1"' \ ' }},' \ ' {{' \ - ' "a_srs":"EPSG:{epsg}",' \ + ' "a_srs":"{srs}",' \ ' "offset_x":"{east}",' \ ' "offset_y":"{north}",' \ ' "offset_z":"0",' \ @@ -225,14 +270,14 @@ class ODM_GeoRef(object): gcp = self.gcps[idx] - kwargs = {'epsg': self.epsg, + kwargs = {'proj': self.projection, 'file': _file, 'x': gcp.x + self.utm_east_offset, 'y': gcp.y + self.utm_north_offset, 'z': gcp.z} latlon = system.run_and_return('echo {x} {y} {z} '.format(**kwargs), - 'gdaltransform -s_srs \"EPSG:{epsg}\" ' + 'gdaltransform -s_srs \"{proj}\" ' '-t_srs \"EPSG:4326\"'.format(**kwargs)).split() # Example: 83d18'16.285"W @@ -297,43 +342,24 @@ class ODM_GeoRef(object): # # write values # metadata.write() - def parse_coordinate_system(self, _file): - """Write attributes to jobOptions from coord file""" - # check for coordinate file existence + def extract_offsets(self, _file): if not io.file_exists(_file): log.ODM_ERROR('Could not find file %s' % _file) return with open(_file) as f: - # extract reference system and utm zone from first line. - # We will assume the following format: - # 'WGS84 UTM 17N' or 'WGS84 UTM 17N \n' - line = f.readline().rstrip() - log.ODM_DEBUG('Line: %s' % line) - ref = line.split(' ') - # match_wgs_utm = re.search('WGS84 UTM (\d{1,2})(N|S)', line, re.I) - if ref[0] == 'WGS84' and ref[1] == 'UTM': # match_wgs_utm: - self.datum = ref[0] - self.utm_pole = ref[2][len(ref[2]) - 1] - self.utm_zone = int(ref[2][:len(ref[2]) - 1]) - # extract east and west offsets from second line. - # We will assume the following format: - # '440143 4588391' - # update EPSG - self.epsg = self.calculate_EPSG(self.utm_zone, self.utm_pole) - # If the first line looks like "EPSG:n" or "epsg:n" - elif ref[0].split(':')[0].lower() == 'epsg': - self.epsg = line.split(':')[1] - else: - log.ODM_ERROR('Could not parse coordinates. Bad CRS supplied: %s' % line) - return - - offsets = f.readline().split(' ') + offsets = f.readlines()[1].split(' ') self.utm_east_offset = int(offsets[0]) self.utm_north_offset = int(offsets[1]) + def create_gcps(self, _file): + if not io.file_exists(_file): + log.ODM_ERROR('Could not find file %s' % _file) + return + + with open(_file) as f: # parse coordinates - lines = f.readlines() + lines = f.readlines()[2:] for l in lines: xyz = l.split(' ') if len(xyz) == 3: @@ -344,6 +370,19 @@ class ODM_GeoRef(object): self.gcps.append(ODM_GCPoint(float(x), float(y), float(z))) # Write to json file + def parse_transformation_matrix(self, _file): + if not io.file_exists(_file): + log.ODM_ERROR('Could not find file %s' % _file) + return + + # Create a nested list for the transformation matrix + with open(_file) as f: + for line in f: + self.transform += [[float(i) for i in line.split()]] + + self.utm_east_offset = self.transform[0][3] + self.utm_north_offset = self.transform[1][3] + class ODM_Tree(object): def __init__(self, root_path, images_path): @@ -374,6 +413,7 @@ class ODM_Tree(object): # benchmarking self.benchmarking = io.join_paths(self.root_path, 'benchmark.txt') + self.dataset_list = io.join_paths(self.root_path, 'img_list.txt') # opensfm self.opensfm_tracks = io.join_paths(self.opensfm, 'tracks.csv') @@ -384,6 +424,7 @@ class ODM_Tree(object): self.opensfm_reconstruction_meshed = io.join_paths(self.opensfm, 'reconstruction.meshed.json') self.opensfm_reconstruction_nvm = io.join_paths(self.opensfm, 'reconstruction.nvm') self.opensfm_model = io.join_paths(self.opensfm, 'depthmaps/merged.ply') + self.opensfm_transformation = io.join_paths(self.opensfm, 'geocoords_transformation.txt') # pmvs self.pmvs_rec_path = io.join_paths(self.pmvs, 'recon0') @@ -410,9 +451,8 @@ class ODM_Tree(object): self.odm_georeferencing_latlon = io.join_paths( self.odm_georeferencing, 'latlon.txt') self.odm_georeferencing_coords = io.join_paths( - self.odm_georeferencing, 'coords.txt') - self.odm_georeferencing_gcp = io.join_paths( - self.odm_georeferencing, 'gcp_list.txt') + self.root_path, 'coords.txt') # Todo put this somewhere better + self.odm_georeferencing_gcp = io.find('gcp_list.txt', self.root_path) self.odm_georeferencing_utm_log = io.join_paths( self.odm_georeferencing, 'odm_georeferencing_utm_log.txt') self.odm_georeferencing_log = 'odm_georeferencing_log.txt' diff --git a/scripts/dataset.py b/scripts/dataset.py index 54285ae1..5f11d43d 100644 --- a/scripts/dataset.py +++ b/scripts/dataset.py @@ -1,8 +1,6 @@ import os import ecto -from functools import partial -from multiprocessing import Pool from opendm import context from opendm import io from opendm import types @@ -24,10 +22,12 @@ class ODMLoadDatasetCell(ecto.Cell): 'images', None) params.declare("force_ccd", 'Override the ccd width information for the ' 'images', None) + params.declare("verbose", 'indicate verbosity', False) + params.declare("proj", 'Geographic projection', None) def declare_io(self, params, inputs, outputs): inputs.declare("tree", "Struct with paths", []) - outputs.declare("photos", "list of ODMPhotos", []) + outputs.declare("reconstruction", "ODMReconstruction", []) def process(self, inputs, outputs): # check if the extension is supported @@ -62,22 +62,48 @@ class ODMLoadDatasetCell(ecto.Cell): if files: # create ODMPhoto list path_files = [io.join_paths(images_dir, f) for f in files] - # photos = Pool().map( - # partial(make_odm_photo, self.params.force_focal, self.params.force_ccd), - # path_files - # ) photos = [] - for files in path_files: - photos += [make_odm_photo(self.params.force_focal, self.params.force_ccd, files)] - - log.ODM_INFO('Found %s usable images' % len(photos)) + with open(tree.dataset_list, 'w') as dataset_list: + for files in path_files: + photos += [make_odm_photo(self.params.force_focal, self.params.force_ccd, files)] + dataset_list.write(photos[-1].filename + '\n') + + log.ODM_INFO('Found %s usable images' % len(photos)) else: log.ODM_ERROR('Not enough supported images in %s' % images_dir) return ecto.QUIT # append photos to cell output - outputs.photos = photos + if not self.params.proj: + if tree.odm_georeferencing_gcp: + outputs.reconstruction = types.ODM_Reconstruction(photos, coords_file=tree.odm_georeferencing_gcp) + else: + verbose = '-verbose' if self.params.verbose else '' + # Generate UTM from images + # odm_georeference definitions + kwargs = { + 'bin': context.odm_modules_path, + 'imgs': tree.dataset_raw, + 'imgs_list': tree.dataset_list, + 'coords': tree.odm_georeferencing_coords, + 'log': tree.odm_georeferencing_utm_log, + 'verbose': verbose + } + + # run UTM extraction binary + extract_utm = system.run_and_return('{bin}/odm_extract_utm -imagesPath {imgs}/ ' + '-imageListFile {imgs_list} -outputCoordFile {coords} {verbose} ' + '-logFile {log}'.format(**kwargs)) + + if extract_utm != '': + log.ODM_WARNING('Could not generate coordinates file. ' + 'Ignore if there is a GCP file. Error: %s' + % extract_utm) + + outputs.reconstruction = types.ODM_Reconstruction(photos, coords_file=tree.odm_georeferencing_coords) + else: + outputs.reconstruction = types.ODM_Reconstruction(photos, projstring=self.params.proj) log.ODM_INFO('Running ODM Load Dataset Cell - Finished') return ecto.OK diff --git a/scripts/mvstex.py b/scripts/mvstex.py index a5a59733..97c8cf21 100644 --- a/scripts/mvstex.py +++ b/scripts/mvstex.py @@ -24,8 +24,6 @@ class ODMMvsTexCell(ecto.Cell): inputs.declare("reconstruction", "Clusters output. list of ODMReconstructions", []) outputs.declare("reconstruction", "Clusters output. list of ODMReconstructions", []) - - def process(self, inputs, outputs): # Benchmarking @@ -34,8 +32,9 @@ class ODMMvsTexCell(ecto.Cell): log.ODM_INFO('Running MVS Texturing Cell') # get inputs - args = self.inputs.args - tree = self.inputs.tree + args = inputs.args + tree = inputs.tree + reconstruction = inputs.reconstruction # define paths and create working directories system.mkdir_p(tree.odm_texturing) @@ -65,8 +64,7 @@ class ODMMvsTexCell(ecto.Cell): if not io.file_exists(odm_textured_model_obj) or rerun_cell: log.ODM_DEBUG('Writing MVS Textured file in: %s' % odm_textured_model_obj) - - + # Format arguments to fit Mvs-Texturing app skipGeometricVisibilityTest = "" skipGlobalSeamLeveling = "" @@ -126,6 +124,8 @@ class ODMMvsTexCell(ecto.Cell): log.ODM_WARNING('Found a valid ODM Texture file in: %s' % odm_textured_model_obj) + outputs.reconstruction = reconstruction + if args.time: system.benchmark(start_time, tree.benchmarking, 'Texturing') diff --git a/scripts/odm_app.py b/scripts/odm_app.py index 1bc5aa24..462c8d5d 100644 --- a/scripts/odm_app.py +++ b/scripts/odm_app.py @@ -7,7 +7,7 @@ from opendm import io from opendm import system from dataset import ODMLoadDatasetCell -from opensfm import ODMOpenSfMCell +from run_opensfm import ODMOpenSfMCell from odm_slam import ODMSlamCell from pmvs import ODMPmvsCell from cmvs import ODMCmvsCell @@ -38,7 +38,9 @@ class ODMApp(ecto.BlackBox): """ cells = {'args': ecto.Constant(value=p.args), 'dataset': ODMLoadDatasetCell(force_focal=p.args.force_focal, - force_ccd=p.args.force_ccd), + force_ccd=p.args.force_ccd, + verbose=p.args.verbose, + proj=p.args.proj), 'opensfm': ODMOpenSfMCell(use_exif_size=False, feature_process_size=p.args.resize_to, feature_min_frames=p.args.min_num_features, @@ -113,7 +115,7 @@ class ODMApp(ecto.BlackBox): # run opensfm with images from load dataset connections += [self.tree[:] >> self.opensfm['tree'], self.args[:] >> self.opensfm['args'], - self.dataset['photos'] >> self.opensfm['photos']] + self.dataset['reconstruction'] >> self.opensfm['reconstruction']] if not p.args.use_pmvs: # create odm mesh from opensfm point cloud @@ -144,7 +146,6 @@ class ODMApp(ecto.BlackBox): # create odm georeference connections += [self.tree[:] >> self.georeferencing['tree'], self.args[:] >> self.georeferencing['args'], - self.dataset['photos'] >> self.georeferencing['photos'], self.texturing['reconstruction'] >> self.georeferencing['reconstruction']] # create odm dem diff --git a/scripts/odm_georeferencing.py b/scripts/odm_georeferencing.py index 0ca6e167..99eec8b9 100644 --- a/scripts/odm_georeferencing.py +++ b/scripts/odm_georeferencing.py @@ -21,64 +21,27 @@ class ODMGeoreferencingCell(ecto.Cell): def declare_io(self, params, inputs, outputs): inputs.declare("tree", "Struct with paths", []) inputs.declare("args", "The application arguments.", {}) - inputs.declare("photos", "list of ODMPhoto's", []) inputs.declare("reconstruction", "list of ODMReconstructions", []) outputs.declare("reconstruction", "list of ODMReconstructions", []) def process(self, inputs, outputs): - # find a file in the root directory - def find(file, dir): - for root, dirs, files in os.walk(dir): - return '/'.join((root, file)) if file in files else None - # Benchmarking start_time = system.now_raw() log.ODM_INFO('Running ODM Georeferencing Cell') # get inputs - args = self.inputs.args - tree = self.inputs.tree - gcpfile = io.join_paths(tree.root_path, self.params.gcp_file) \ - if self.params.gcp_file else find('gcp_list.txt', tree.root_path) + args = inputs.args + tree = inputs.tree + reconstruction = inputs.reconstruction + gcpfile = tree.odm_georeferencing_gcp geocreated = True verbose = '-verbose' if self.params.verbose else '' # define paths and create working directories system.mkdir_p(tree.odm_georeferencing) if args.use_25dmesh: system.mkdir_p(tree.odm_25dgeoreferencing) - - # in case a gcp file it's not provided, let's try to generate it using - # images metadata. Internally calls jhead. - log.ODM_DEBUG(self.params.gcp_file) - if not self.params.gcp_file: # and \ - # not io.file_exists(tree.odm_georeferencing_coords): - - log.ODM_WARNING('No coordinates file. ' - 'Generating coordinates file: %s' - % tree.odm_georeferencing_coords) - - # odm_georeference definitions - kwargs = { - 'bin': context.odm_modules_path, - 'imgs': tree.dataset_raw, - 'imgs_list': tree.opensfm_bundle_list, - 'coords': tree.odm_georeferencing_coords, - 'log': tree.odm_georeferencing_utm_log, - 'verbose': verbose - } - - # run UTM extraction binary - extract_utm = system.run_and_return('{bin}/odm_extract_utm -imagesPath {imgs}/ ' - '-imageListFile {imgs_list} -outputCoordFile {coords} {verbose} ' - '-logFile {log}'.format(**kwargs)) - - if extract_utm != '': - log.ODM_WARNING('Could not generate coordinates file. ' - 'Ignore if there is a GCP file. Error: %s' - % extract_utm) - # check if we rerun cell or not rerun_cell = (args.rerun is not None and @@ -116,6 +79,7 @@ class ODMGeoreferencingCell(ecto.Cell): 'imgs_list': tree.opensfm_bundle_list, 'model': r['model'], 'log': odm_georeferencing_log, + 'input_trans_file': tree.opensfm_transformation, 'transform_file': odm_georeferencing_transform_file, 'coords': tree.odm_georeferencing_coords, 'pc_geo': odm_georeferencing_model_ply_geo, @@ -132,17 +96,23 @@ class ODMGeoreferencingCell(ecto.Cell): # Check to see if the GCP file exists - if not self.params.use_exif and (self.params.gcp_file or find('gcp_list.txt', tree.root_path)): - log.ODM_INFO('Found %s' % gcpfile) - try: - system.run('{bin}/odm_georef -bundleFile {bundle} -imagesPath {imgs} -imagesListPath {imgs_list} ' - '-inputFile {model} -outputFile {model_geo} ' - '-inputPointCloudFile {pc} -outputPointCloudFile {pc_geo} {verbose} ' - '-logFile {log} -outputTransformFile {transform_file} -georefFileOutputPath {geo_sys} -gcpFile {gcp} ' - '-outputCoordFile {coords}'.format(**kwargs)) - except Exception: - log.ODM_EXCEPTION('Georeferencing failed. ') - return ecto.QUIT + #if not self.params.use_exif and (self.params.gcp_file or tree.odm_georeferencing_gcp): + # log.ODM_INFO('Found %s' % gcpfile) + # try: + # system.run('{bin}/odm_georef -bundleFile {bundle} -imagesPath {imgs} -imagesListPath {imgs_list} ' + # '-inputFile {model} -outputFile {model_geo} ' + # '-inputPointCloudFile {pc} -outputPointCloudFile {pc_geo} {verbose} ' + # '-logFile {log} -outputTransformFile {transform_file} -georefFileOutputPath {geo_sys} -gcpFile {gcp} ' + # '-outputCoordFile {coords}'.format(**kwargs)) + # except Exception: + # log.ODM_EXCEPTION('Georeferencing failed. ') + # return ecto.QUIT + if io.file_exists(tree.opensfm_transformation): + log.ODM_INFO('Running georeferencing with OpenSfM transformation matrix') + system.run('{bin}/odm_georef -bundleFile {bundle} -inputTransformFile {input_trans_file} ' + '-inputFile {model} -outputFile {model_geo} ' + '-inputPointCloudFile {pc} -outputPointCloudFile {pc_geo} {verbose} ' + '-logFile {log} -outputTransformFile {transform_file} -georefFileOutputPath {geo_sys}'.format(**kwargs)) elif io.file_exists(tree.odm_georeferencing_coords): log.ODM_INFO('Running georeferencing with generated coords file.') system.run('{bin}/odm_georef -bundleFile {bundle} -inputCoordFile {coords} ' @@ -153,22 +123,25 @@ class ODMGeoreferencingCell(ecto.Cell): log.ODM_WARNING('Georeferencing failed. Make sure your ' 'photos have geotags in the EXIF or you have ' 'provided a GCP file. ') - geocreated = False # skip the rest of the georeferencing + geocreated = False # skip the rest of the georeferencing odm_georeferencing_model_ply_geo = os.path.join(tree.odm_georeferencing, tree.odm_georeferencing_model_ply_geo) if geocreated: # update images metadata - geo_ref = types.ODM_GeoRef() - geo_ref.parse_coordinate_system(tree.odm_georeferencing_coords) + geo_ref = reconstruction.georef + geo_ref.extract_offsets(tree.odm_georeferencing_coords) + geo_ref.parse_transformation_matrix(tree.opensfm_transformation) - for idx, photo in enumerate(self.inputs.photos): - geo_ref.utm_to_latlon(tree.odm_georeferencing_latlon, photo, idx) + # for idx, photo in enumerate(reconstruction.photos): + # geo_ref.utm_to_latlon(tree.odm_georeferencing_latlon, photo, idx) # convert ply model to LAS reference system geo_ref.convert_to_las(odm_georeferencing_model_ply_geo, tree.odm_georeferencing_model_las, tree.odm_georeferencing_las_json) + reconstruction.georef = geo_ref + # XYZ point cloud output log.ODM_INFO("Creating geo-referenced CSV file (XYZ format)") with open(tree.odm_georeferencing_xyz_file, "wb") as csvfile: @@ -190,6 +163,8 @@ class ODMGeoreferencingCell(ecto.Cell): log.ODM_WARNING('Found a valid georeferenced model in: %s' % odm_georeferencing_model_ply_geo) + outputs.reconstruction = reconstruction + if args.time: system.benchmark(start_time, tree.benchmarking, 'Georeferencing') diff --git a/scripts/odm_meshing.py b/scripts/odm_meshing.py index 8a05552c..520ad410 100644 --- a/scripts/odm_meshing.py +++ b/scripts/odm_meshing.py @@ -42,8 +42,9 @@ class ODMeshingCell(ecto.Cell): log.ODM_INFO('Running ODM Meshing Cell') # get inputs - args = self.inputs.args - tree = self.inputs.tree + args = inputs.args + tree = inputs.tree + reconstruction = inputs.reconstruction verbose = '-verbose' if self.params.verbose else '' # define paths and create working directories @@ -109,6 +110,8 @@ class ODMeshingCell(ecto.Cell): log.ODM_WARNING('Found a valid ODM 2.5D Mesh file in: %s' % tree.odm_25dmesh) + outputs.reconstruction = reconstruction + if args.time: system.benchmark(start_time, tree.benchmarking, 'Meshing') diff --git a/scripts/odm_orthophoto.py b/scripts/odm_orthophoto.py index 4174af42..825a9a27 100644 --- a/scripts/odm_orthophoto.py +++ b/scripts/odm_orthophoto.py @@ -32,6 +32,7 @@ class ODMOrthoPhotoCell(ecto.Cell): # get inputs args = self.inputs.args tree = self.inputs.tree + reconstruction = inputs.reconstruction verbose = '-verbose' if self.params.verbose else '' # define paths and create working directories @@ -57,7 +58,7 @@ class ODMOrthoPhotoCell(ecto.Cell): } # Have geo coordinates? - if io.file_exists(tree.odm_georeferencing_coords): + if reconstruction.georef: # io.file_exists(tree.odm_georeferencing_coords): if args.use_25dmesh: kwargs['model_geo'] = os.path.join(tree.odm_25dtexturing, tree.odm_georeferencing_model_obj_geo) else: @@ -79,11 +80,11 @@ class ODMOrthoPhotoCell(ecto.Cell): else: # Create georeferenced GeoTiff geotiffcreated = False - georef = types.ODM_GeoRef() + georef = reconstruction.georef # creates the coord refs # TODO I don't want to have to do this twice- after odm_georef - georef.parse_coordinate_system(tree.odm_georeferencing_coords) + # georef.parse_coordinate_system(tree.odm_georeferencing_coords) - if georef.epsg and georef.utm_east_offset and georef.utm_north_offset: + if georef.projection and georef.utm_east_offset and georef.utm_north_offset: ulx = uly = lrx = lry = 0.0 with open(tree.odm_orthophoto_corners) as f: for lineNumber, line in enumerate(f): @@ -109,8 +110,7 @@ class ODMOrthoPhotoCell(ecto.Cell): 'compress': self.params.compress, 'predictor': '-co PREDICTOR=2 ' if self.params.compress in ['LZW', 'DEFLATE'] else '', - 'epsg': georef.epsg, - 't_srs': self.params.t_srs or "EPSG:{0}".format(georef.epsg), + 'proj': georef.projection.srs, 'bigtiff': self.params.bigtiff, 'png': tree.odm_orthophoto_file, 'tiff': tree.odm_orthophoto_tif, @@ -125,7 +125,7 @@ class ODMOrthoPhotoCell(ecto.Cell): '-co BLOCKXSIZE=512 ' '-co BLOCKYSIZE=512 ' '-co NUM_THREADS=ALL_CPUS ' - '-a_srs \"EPSG:{epsg}\" ' + '-a_srs \"{proj}\" ' '{png} {tiff} > {log}'.format(**kwargs)) if self.params.build_overviews: diff --git a/scripts/opensfm.py b/scripts/run_opensfm.py similarity index 92% rename from scripts/opensfm.py rename to scripts/run_opensfm.py index 76aa28b3..090e06da 100644 --- a/scripts/opensfm.py +++ b/scripts/run_opensfm.py @@ -20,7 +20,7 @@ class ODMOpenSfMCell(ecto.Cell): def declare_io(self, params, inputs, outputs): inputs.declare("tree", "Struct with paths", []) inputs.declare("args", "The application arguments.", {}) - inputs.declare("photos", "list of ODMPhoto's", []) + inputs.declare("reconstruction", "ODMReconstruction", []) outputs.declare("reconstruction", "list of ODMReconstructions", []) def process(self, inputs, outputs): @@ -31,9 +31,10 @@ class ODMOpenSfMCell(ecto.Cell): log.ODM_INFO('Running ODM OpenSfM Cell') # get inputs - tree = self.inputs.tree - args = self.inputs.args - photos = self.inputs.photos + tree = inputs.tree + args = inputs.args + reconstruction = inputs.reconstruction + photos = reconstruction.photos if not photos: log.ODM_ERROR('Not enough photos in photos array to start OpenSfM') @@ -73,6 +74,7 @@ class ODMOpenSfMCell(ecto.Cell): "feature_min_frames: %s" % self.params.feature_min_frames, "processes: %s" % self.params.processes, "matching_gps_neighbors: %s" % self.params.matching_gps_neighbors, + "depthmap_resolution: 640", "optimize_camera_parameters: %s" % ('no' if self.params.fixed_camera_params else 'yes') ] @@ -90,7 +92,12 @@ class ODMOpenSfMCell(ecto.Cell): if args.matcher_distance > 0: config.append("matching_gps_distance: %s" % self.params.matching_gps_distance) + if tree.odm_georeferencing_gcp: + config.append("bundle_use_gcp: yes") + io.copy(tree.odm_georeferencing_gcp, tree.opensfm) + # write config file + log.ODM_DEBUG(config) config_filename = io.join_paths(tree.opensfm, 'config.yaml') with open(config_filename, 'w') as fout: fout.write("\n".join(config)) @@ -165,6 +172,11 @@ class ODMOpenSfMCell(ecto.Cell): else: log.ODM_WARNING('Found a valid CMVS file in: %s' % tree.pmvs_visdat) + system.run('PYTHONPATH=%s %s/bin/opensfm export_geocoords %s --transformation --proj \'%s\'' % + (context.pyopencv_path, context.opensfm_path, tree.opensfm, reconstruction.georef.projection.srs)) + + outputs.reconstruction = reconstruction + if args.time: system.benchmark(start_time, tree.benchmarking, 'OpenSfM')