kopia lustrzana https://github.com/OpenDroneMap/ODM
Refactored odm_orthophoto to use GDAL for TIFF writing
Piero Toffanin
rodzic
f385d39a9a
commit
f484cc6771
|
|
@ -11,6 +11,8 @@ add_definitions(-Wall -Wextra)
|
|||
# Find pcl at the location specified by PCL_DIR
|
||||
find_package(VTK 6.0 REQUIRED)
|
||||
find_package(PCL 1.8 HINTS "${PCL_DIR}/share/pcl-1.8" REQUIRED)
|
||||
find_package(GDAL REQUIRED)
|
||||
include_directories(${GDAL_INCLUDE_DIR})
|
||||
|
||||
# Find OpenCV at the default location
|
||||
find_package(OpenCV HINTS "${OPENCV_DIR}" REQUIRED)
|
||||
|
|
@ -31,4 +33,4 @@ aux_source_directory("./src" SRC_LIST)
|
|||
|
||||
# Add exectuteable
|
||||
add_executable(${PROJECT_NAME} ${SRC_LIST})
|
||||
target_link_libraries(odm_orthophoto ${PCL_COMMON_LIBRARIES} ${PCL_IO_LIBRARIES} ${PCL_SURFACE_LIBRARIES} ${OpenCV_LIBS})
|
||||
target_link_libraries(odm_orthophoto ${PCL_COMMON_LIBRARIES} ${PCL_IO_LIBRARIES} ${PCL_SURFACE_LIBRARIES} ${OpenCV_LIBS} ${GDAL_LIBRARY})
|
||||
|
|
|
|||
|
|
@ -306,6 +306,50 @@ void OdmOrthoPhoto::printHelp()
|
|||
log_.setIsPrintingInCout(false);
|
||||
}
|
||||
|
||||
void OdmOrthoPhoto::saveTIFF(const std::string &filename, GDALDataType dataType){
|
||||
GDALAllRegister();
|
||||
GDALDriverH hDriver = GDALGetDriverByName( "GTiff" );
|
||||
if (!hDriver){
|
||||
std::cerr << "Cannot initialize GeoTIFF driver. Check your GDAL installation." << std::endl;
|
||||
exit(1);
|
||||
}
|
||||
char **papszOptions = NULL;
|
||||
GDALDatasetH hDstDS = GDALCreate( hDriver, filename.c_str(), width, height, static_cast<int>(bands.size()), dataType, papszOptions );
|
||||
GDALRasterBandH hBand;
|
||||
|
||||
for (size_t i = 0; i < bands.size(); i++){
|
||||
hBand = GDALGetRasterBand( hDstDS, static_cast<int>(i) + 1 );
|
||||
if (GDALRasterIO( hBand, GF_Write, 0, 0, width, height,
|
||||
bands[i], width, height, dataType, 0, 0 ) != CE_None){
|
||||
std::cerr << "Cannot write TIFF to " << filename << std::endl;
|
||||
exit(1);
|
||||
}
|
||||
// for (int j = 0; j < height; j++){
|
||||
// GDALRasterIO( hBand, GF_Write, 0, j, width, 1,
|
||||
// bands[i][j], width, 1, dataType, 0, 0 );
|
||||
// }
|
||||
}
|
||||
GDALClose( hDstDS );
|
||||
}
|
||||
|
||||
template <typename T>
|
||||
inline T maxRange(){
|
||||
return static_cast<T>(pow(2, sizeof(T) * 8) - 1);
|
||||
}
|
||||
|
||||
template <typename T>
|
||||
void OdmOrthoPhoto::initBands(const cv::Mat &texture){
|
||||
// Channels + alpha
|
||||
for (int i = 0; i < texture.channels() + 1; i++){
|
||||
size_t pixelCount = static_cast<size_t>(width * height);
|
||||
T *arr = new T[pixelCount];
|
||||
for (size_t j = 0; j < pixelCount; j++){
|
||||
arr[j] = i == 4 ? 0 : maxRange<T>(); // Set alpha at 0
|
||||
}
|
||||
bands.push_back(static_cast<void *>(arr));
|
||||
}
|
||||
}
|
||||
|
||||
void OdmOrthoPhoto::createOrthoPhoto()
|
||||
{
|
||||
if(inputFile_.empty())
|
||||
|
|
@ -333,7 +377,7 @@ void OdmOrthoPhoto::createOrthoPhoto()
|
|||
log_ << "\tSpecified -inputGeoRefFile, but no boundary points. The georeference system will be ignored.\n";
|
||||
}
|
||||
|
||||
log_ << "Reading mesh file...\n";
|
||||
log_ << "Reading mesh file... " << inputFile_ << "\n";
|
||||
// The textured mesh.
|
||||
pcl::TextureMesh mesh;
|
||||
loadObjFile(inputFile_, mesh);
|
||||
|
|
@ -375,24 +419,24 @@ void OdmOrthoPhoto::createOrthoPhoto()
|
|||
log_ << "Ortho photo area : " << xDiff*yDiff << "m2\n";
|
||||
|
||||
// The resolution necessary to fit the area with the given resolution.
|
||||
int rowRes = static_cast<int>(std::ceil(resolution_*yDiff));
|
||||
int colRes = static_cast<int>(std::ceil(resolution_*xDiff));
|
||||
log_ << "Ortho photo resolution, width x height : " << colRes << "x" << rowRes << '\n';
|
||||
height = static_cast<int>(std::ceil(resolution_*yDiff));
|
||||
width = static_cast<int>(std::ceil(resolution_*xDiff));
|
||||
log_ << "Ortho photo resolution, width x height : " << width << "x" << height << '\n';
|
||||
|
||||
// Check size of photo.
|
||||
if(0 >= rowRes*colRes)
|
||||
if(0 >= height*width)
|
||||
{
|
||||
if(0 >= rowRes)
|
||||
if(0 >= height)
|
||||
{
|
||||
log_ << "Warning: ortho photo has zero area, height = " << rowRes << ". Forcing height = 1.\n";
|
||||
rowRes = 1;
|
||||
log_ << "Warning: ortho photo has zero area, height = " << height << ". Forcing height = 1.\n";
|
||||
height = 1;
|
||||
}
|
||||
if(0 >= colRes)
|
||||
if(0 >= width)
|
||||
{
|
||||
log_ << "Warning: ortho photo has zero area, width = " << colRes << ". Forcing width = 1.\n";
|
||||
colRes = 1;
|
||||
log_ << "Warning: ortho photo has zero area, width = " << width << ". Forcing width = 1.\n";
|
||||
width = 1;
|
||||
}
|
||||
log_ << "New ortho photo resolution, width x height : " << colRes << "x" << rowRes << '\n';
|
||||
log_ << "New ortho photo resolution, width x height : " << width << "x" << height << '\n';
|
||||
}
|
||||
|
||||
// Contains the vertices of the mesh.
|
||||
|
|
@ -493,20 +537,22 @@ void OdmOrthoPhoto::createOrthoPhoto()
|
|||
if (textureDepth == -1){
|
||||
try{
|
||||
textureDepth = texture.depth();
|
||||
log_ << "Texture depth is " << textureDepth << "bit\n";
|
||||
log_ << "Texture channels: " << texture.channels() << "\n";
|
||||
|
||||
if (textureDepth == CV_8U){
|
||||
photo_ = cv::Mat::zeros(rowRes, colRes, CV_8UC4) + cv::Scalar(255, 255, 255, 0);
|
||||
log_ << "Texture depth: 8bit\n";
|
||||
initBands<uint8_t>(texture);
|
||||
}else if (textureDepth == CV_16U){
|
||||
photo_ = cv::Mat::zeros(rowRes, colRes, CV_16UC4) + cv::Scalar(65535, 65535, 65535, 0);
|
||||
log_ << "Texture depth: 16bit\n";
|
||||
initBands<uint16_t>(texture);
|
||||
}else{
|
||||
std::cerr << "Unsupported bit depth value: " << textureDepth;
|
||||
exit(1);
|
||||
}
|
||||
|
||||
depth_ = cv::Mat::zeros(rowRes, colRes, CV_32F) - std::numeric_limits<float>::infinity();
|
||||
depth_ = cv::Mat::zeros(height, width, CV_32F) - std::numeric_limits<float>::infinity();
|
||||
}catch(const cv::Exception &e){
|
||||
std::cerr << "Couldn't allocate enough memory to render the orthophoto (" << colRes << "x" << rowRes << " cells = " << ((long long)colRes * (long long)rowRes * 4) << " bytes). Try to increase the --orthophoto-resolution parameter to a larger integer or add more RAM.\n";
|
||||
std::cerr << "Couldn't allocate enough memory to render the orthophoto (" << width << "x" << height << " cells = " << ((long long)width * (long long)height * 4) << " bytes). Try to increase the --orthophoto-resolution parameter to a larger integer or add more RAM.\n";
|
||||
exit(1);
|
||||
}
|
||||
}
|
||||
|
|
@ -531,9 +577,9 @@ void OdmOrthoPhoto::createOrthoPhoto()
|
|||
|
||||
// ... and draw it into the ortho photo.
|
||||
if (textureDepth == CV_8U){
|
||||
drawTexturedTriangle<unsigned char>(texture, polygon, meshCloud, uvs, faceIndex+faceOff);
|
||||
drawTexturedTriangle<uint8_t>(texture, polygon, meshCloud, uvs, faceIndex+faceOff);
|
||||
}else if (textureDepth == CV_16U){
|
||||
drawTexturedTriangle<unsigned short>(texture, polygon, meshCloud, uvs, faceIndex+faceOff);
|
||||
drawTexturedTriangle<uint16_t>(texture, polygon, meshCloud, uvs, faceIndex+faceOff);
|
||||
}
|
||||
}
|
||||
faceOff += faces.size();
|
||||
|
|
@ -543,7 +589,17 @@ void OdmOrthoPhoto::createOrthoPhoto()
|
|||
|
||||
log_ << '\n';
|
||||
log_ << "Writing ortho photo to " << outputFile_ << "\n";
|
||||
cv::imwrite(outputFile_, photo_);
|
||||
|
||||
if (textureDepth == CV_8U){
|
||||
saveTIFF(outputFile_, GDT_Byte);
|
||||
}else if (textureDepth == CV_16U){
|
||||
saveTIFF(outputFile_, GDT_UInt16);
|
||||
}else{
|
||||
std::cerr << "Unsupported bit depth value: " << textureDepth;
|
||||
exit(1);
|
||||
}
|
||||
|
||||
// cv::imwrite(outputFile_, photo_);
|
||||
|
||||
if (!outputCornerFile_.empty())
|
||||
{
|
||||
|
|
@ -801,7 +857,7 @@ void OdmOrthoPhoto::drawTexturedTriangle(const cv::Mat &texture, const pcl::Vert
|
|||
|
||||
// Check bounding box overlap.
|
||||
int xMin = static_cast<int>(std::min(std::min(v1x, v2x), v3x));
|
||||
if(xMin > photo_.cols)
|
||||
if(xMin > width)
|
||||
{
|
||||
return; // Completely outside to the right.
|
||||
}
|
||||
|
|
@ -811,7 +867,7 @@ void OdmOrthoPhoto::drawTexturedTriangle(const cv::Mat &texture, const pcl::Vert
|
|||
return; // Completely outside to the left.
|
||||
}
|
||||
int yMin = static_cast<int>(std::min(std::min(v1y, v2y), v3y));
|
||||
if(yMin > photo_.rows)
|
||||
if(yMin > height)
|
||||
{
|
||||
return; // Completely outside to the top.
|
||||
}
|
||||
|
|
@ -915,7 +971,7 @@ void OdmOrthoPhoto::drawTexturedTriangle(const cv::Mat &texture, const pcl::Vert
|
|||
// The first pixel row for the bottom part of the triangle.
|
||||
int rqStart = std::max(static_cast<int>(std::floor(topR+0.5f)), 0);
|
||||
// The last pixel row for the top part of the triangle.
|
||||
int rqEnd = std::min(static_cast<int>(std::floor(midR+0.5f)), photo_.rows);
|
||||
int rqEnd = std::min(static_cast<int>(std::floor(midR+0.5f)), height);
|
||||
|
||||
// Traverse along row from top to middle.
|
||||
for(int rq = rqStart; rq < rqEnd; ++rq)
|
||||
|
|
@ -927,7 +983,7 @@ void OdmOrthoPhoto::drawTexturedTriangle(const cv::Mat &texture, const pcl::Vert
|
|||
// The first pixel column for the current row.
|
||||
int cqStart = std::max(static_cast<int>(std::floor(0.5f+std::min(ctm, ctb))), 0);
|
||||
// The last pixel column for the current row.
|
||||
int cqEnd = std::min(static_cast<int>(std::floor(0.5f+std::max(ctm, ctb))), photo_.cols);
|
||||
int cqEnd = std::min(static_cast<int>(std::floor(0.5f+std::max(ctm, ctb))), width);
|
||||
|
||||
for(int cq = cqStart; cq < cqEnd; ++cq)
|
||||
{
|
||||
|
|
@ -973,7 +1029,7 @@ void OdmOrthoPhoto::drawTexturedTriangle(const cv::Mat &texture, const pcl::Vert
|
|||
// The first pixel row for the bottom part of the triangle.
|
||||
int rqStart = std::max(static_cast<int>(std::floor(midR+0.5f)), 0);
|
||||
// The last pixel row for the bottom part of the triangle.
|
||||
int rqEnd = std::min(static_cast<int>(std::floor(botR+0.5f)), photo_.rows);
|
||||
int rqEnd = std::min(static_cast<int>(std::floor(botR+0.5f)), height);
|
||||
|
||||
// Traverse along row from middle to bottom.
|
||||
for(int rq = rqStart; rq < rqEnd; ++rq)
|
||||
|
|
@ -985,7 +1041,7 @@ void OdmOrthoPhoto::drawTexturedTriangle(const cv::Mat &texture, const pcl::Vert
|
|||
// The first pixel column for the current row.
|
||||
int cqStart = std::max(static_cast<int>(std::floor(0.5f+std::min(cmb, ctb))), 0);
|
||||
// The last pixel column for the current row.
|
||||
int cqEnd = std::min(static_cast<int>(std::floor(0.5f+std::max(cmb, ctb))), photo_.cols);
|
||||
int cqEnd = std::min(static_cast<int>(std::floor(0.5f+std::max(cmb, ctb))), width);
|
||||
|
||||
for(int cq = cqStart; cq < cqEnd; ++cq)
|
||||
{
|
||||
|
|
@ -1071,10 +1127,11 @@ void OdmOrthoPhoto::renderPixel(int row, int col, float s, float t, const cv::Ma
|
|||
b += static_cast<float>(bl[0]) * dr * dt;
|
||||
b += static_cast<float>(br[0]) * dl * dt;
|
||||
|
||||
photo_.at<cv::Vec<T, 4> >(row,col) = cv::Vec<T, 4>(static_cast<T>(b),
|
||||
static_cast<T>(g),
|
||||
static_cast<T>(r),
|
||||
static_cast<T>(255)); // Alpha should always be in the 255 range
|
||||
size_t idx = static_cast<size_t>(row * width + col);
|
||||
static_cast<T *>(bands[0])[idx] = static_cast<T>(r);
|
||||
static_cast<T *>(bands[1])[idx] = static_cast<T>(g);
|
||||
static_cast<T *>(bands[2])[idx] = static_cast<T>(b);
|
||||
static_cast<T *>(bands[3])[idx] = static_cast<T>(255); // Alpha should always be in the 255 range
|
||||
}
|
||||
|
||||
void OdmOrthoPhoto::getBarycentricCoordinates(pcl::PointXYZ v1, pcl::PointXYZ v2, pcl::PointXYZ v3, float x, float y, float &l1, float &l2, float &l3) const
|
||||
|
|
|
|||
|
|
@ -20,6 +20,10 @@
|
|||
// OpenCV
|
||||
#include <opencv2/core/core.hpp>
|
||||
|
||||
// GDAL
|
||||
#include "gdal_priv.h"
|
||||
#include "cpl_conv.h" // for CPLMalloc()
|
||||
|
||||
// Logger
|
||||
#include "Logger.hpp"
|
||||
|
||||
|
|
@ -76,6 +80,7 @@ public:
|
|||
int run(int argc, char* argv[]);
|
||||
|
||||
private:
|
||||
int width, height;
|
||||
|
||||
/*!
|
||||
* \brief parseArguments Parses command line arguments.
|
||||
|
|
@ -123,6 +128,10 @@ private:
|
|||
* @return
|
||||
*/
|
||||
Eigen::Transform<float, 3, Eigen::Affine> readTransform(std::string transformFile_) const;
|
||||
|
||||
template <typename T>
|
||||
void initBands(const cv::Mat &texture);
|
||||
void saveTIFF(const std::string &filename, GDALDataType dataType);
|
||||
|
||||
/*!
|
||||
* \brief Renders a triangle into the ortho photo.
|
||||
|
|
@ -222,7 +231,8 @@ private:
|
|||
Eigen::Vector2f boundaryPoint3_; /**< The third boundary point for the ortho photo, in local coordinates. */
|
||||
Eigen::Vector2f boundaryPoint4_; /**< The fourth boundary point for the ortho photo, in local coordinates. */
|
||||
|
||||
cv::Mat photo_; /**< The ortho photo as an OpenCV matrix, CV_8UC3. */
|
||||
std::vector<void *> bands;
|
||||
|
||||
cv::Mat depth_; /**< The depth of the ortho photo as an OpenCV matrix, CV_32F. */
|
||||
|
||||
bool multiMaterial_; /**< True if the mesh has multiple materials. **/
|
||||
|
|
|
|||
|
|
@ -82,11 +82,15 @@ def compute_mask_raster(input_raster, vector_mask, output_raster, blend_distance
|
|||
out_image, out_transform = mask(rast, shapes, nodata=0)
|
||||
|
||||
if blend_distance > 0:
|
||||
alpha_band = out_image[3]
|
||||
dist_t = ndimage.distance_transform_edt(alpha_band)
|
||||
dist_t[dist_t <= blend_distance] /= blend_distance
|
||||
dist_t[dist_t > blend_distance] = 1
|
||||
np.multiply(alpha_band, dist_t, out=alpha_band, casting="unsafe")
|
||||
if out_image.shape[0] >= 4:
|
||||
# rast_mask = rast.dataset_mask()
|
||||
rast_mask = out_image[-1]
|
||||
dist_t = ndimage.distance_transform_edt(rast_mask)
|
||||
dist_t[dist_t <= blend_distance] /= blend_distance
|
||||
dist_t[dist_t > blend_distance] = 1
|
||||
np.multiply(rast_mask, dist_t, out=rast_mask, casting="unsafe")
|
||||
else:
|
||||
log.ODM_WARNING("%s does not have an alpha band, cannot blend cutline!" % input_raster)
|
||||
|
||||
with rasterio.open(output_raster, 'w', **rast.profile) as dst:
|
||||
dst.write(out_image)
|
||||
|
|
|
|||
|
|
@ -93,6 +93,54 @@ class ODM_Reconstruction(object):
|
|||
self.photos = photos
|
||||
self.georef = None
|
||||
self.gcp = None
|
||||
self.multi_camera = self.detect_multi_camera()
|
||||
|
||||
def detect_multi_camera(self):
|
||||
"""
|
||||
Looks at the reconstruction photos and determines if this
|
||||
is a single or multi-camera setup.
|
||||
"""
|
||||
supported_ext_re = r"\.(" + "|".join([e[1:] for e in context.supported_extensions]) + ")$"
|
||||
|
||||
# Match filename_1.tif, filename_2.tif, filename_3.tif, ...
|
||||
#
|
||||
multi_camera_patterns = {
|
||||
'MicaSense RedEdge-M': r'^IMG_\d+_(?P<band>\d{1})',
|
||||
'Parrot Sequoia': r'^IMG_\d+_\d+_\d+_(?P<band>[A-Z]{3})',
|
||||
}
|
||||
|
||||
for cam, regex in multi_camera_patterns.items():
|
||||
pattern = re.compile(regex + supported_ext_re, re.IGNORECASE)
|
||||
mc = {}
|
||||
|
||||
for p in self.photos:
|
||||
matches = re.match(pattern, p.filename)
|
||||
if matches:
|
||||
band = matches.group("band")
|
||||
if not band in mc:
|
||||
mc[band] = []
|
||||
mc[band].append(p)
|
||||
|
||||
# We support between 2 and 6 bands
|
||||
# If we matched more or less bands, we probably just
|
||||
# found filename patterns that do not match a multi-camera setup
|
||||
bands_count = len(mc)
|
||||
if bands_count >= 2 and bands_count <= 6:
|
||||
|
||||
# Validate that all bands have the same number of images,
|
||||
# otherwise this is not a multi-camera setup
|
||||
img_per_band = len(mc[band])
|
||||
valid = True
|
||||
for band in mc:
|
||||
if len(mc[band]) != img_per_band:
|
||||
log.ODM_WARNING("This might be a multi-camera setup, but band \"%s\" (identified from \"%s\") has only %s images (instead of %s), perhaps images are missing or are corrupted." % (band, mc[band][0].filename, len(mc[band]), img_per_band))
|
||||
valid = False
|
||||
break
|
||||
|
||||
if valid:
|
||||
return mc
|
||||
|
||||
return None
|
||||
|
||||
def is_georeferenced(self):
|
||||
return self.georef is not None
|
||||
|
|
|
|||
|
|
@ -97,7 +97,7 @@ class ODMOrthoPhotoStage(types.ODM_Stage):
|
|||
|
||||
system.run('gdal_translate -a_ullr {ulx} {uly} {lrx} {lry} '
|
||||
'{vars} '
|
||||
'-b 1 -b 2 -b 3 -mask 4 '
|
||||
# '-b 1 -b 2 -b 3 -mask 4 '
|
||||
'-a_srs \"{proj}\" '
|
||||
'--config GDAL_CACHEMAX {max_memory}% '
|
||||
'--config GDAL_TIFF_INTERNAL_MASK YES '
|
||||
|
|
|
|||
|
|
@ -35,10 +35,11 @@ class TestRemote(unittest.TestCase):
|
|||
self.queue_num = queue_num
|
||||
self.uuid = 'xxxxx-xxxxx-xxxxx-xxxxx-xxxx' + str(queue_num)
|
||||
|
||||
def info(self):
|
||||
def info(self, with_output=None):
|
||||
class StatusMock:
|
||||
status = TaskStatus.RUNNING if self.running else TaskStatus.QUEUED
|
||||
processing_time = 1
|
||||
output = "test output"
|
||||
return StatusMock()
|
||||
|
||||
def remove(self):
|
||||
|
|
|
|||
Ładowanie…
Reference in New Issue