mirror
/
OpenDroneMap-ODM
kopia lustrzana https://github.com/OpenDroneMap/ODM
1
0
Forkuj 0
Kod Zgłoszenia Projekty Wydania Wiki Aktywność

Merge pull request #904 from pierotofy/dem2mesh 

Dem2mesh
pull/912/head
Piero Toffanin 2018-11-01 10:47:18 -05:00 zatwierdzone przez GitHub
rodzic 002a6f0b1e e2f4a78ddd
commit ee5353350c
Nie znaleziono w bazie danych klucza dla tego podpisu
ID klucza GPG: 4AEE18F83AFDEB23
14 zmienionych plików z 101 dodań i 770 usunięć
Pokaż statystyki Ściągnij plik aktualizacji Ściągnij plik porównania Expand all files Collapse all files

20
SuperBuild/CMakeLists.txt
Wyświetl plik

@ -163,3 +163,23 @@ externalproject_add(poissonrecon
BUILD_COMMAND make poissonrecon
INSTALL_COMMAND ""
)
externalproject_add(dem2mesh
GIT_REPOSITORY https://github.com/OpenDroneMap/dem2mesh.git
GIT_TAG master
SOURCE_DIR ${SB_SOURCE_DIR}/dem2mesh
UPDATE_COMMAND ""
BUILD_IN_SOURCE 1
BUILD_COMMAND make
INSTALL_COMMAND ""
)
externalproject_add(dem2points
GIT_REPOSITORY https://github.com/OpenDroneMap/dem2points.git
GIT_TAG master
SOURCE_DIR ${SB_SOURCE_DIR}/dem2points
UPDATE_COMMAND ""
BUILD_IN_SOURCE 1
BUILD_COMMAND make
INSTALL_COMMAND ""
)

2
VERSION
Wyświetl plik

@ -1 +1 @@
0.4.0
0.4.1

1
modules/CMakeLists.txt
Wyświetl plik

@ -8,7 +8,6 @@ add_subdirectory(odm_extract_utm)
add_subdirectory(odm_georef)
add_subdirectory(odm_orthophoto)
add_subdirectory(odm_cleanmesh)
add_subdirectory(odm_dem2points)
if (ODM_BUILD_SLAM)
add_subdirectory(odm_slam)
endif ()

19
modules/odm_dem2points/CMakeLists.txt
Wyświetl plik

@ -1,19 +0,0 @@
project(odm_dem2points)
cmake_minimum_required(VERSION 2.8)
set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${CMAKE_SOURCE_DIR})
set (CMAKE_CXX_STANDARD 11)
find_package(GDAL REQUIRED)
include_directories(${GDAL_INCLUDE_DIR})
# Add compiler options.
add_definitions(-Wall -Wextra)
# Add source directory
aux_source_directory("./src" SRC_LIST)
# Add exectuteable
add_executable(${PROJECT_NAME} ${SRC_LIST})
target_link_libraries(${PROJECT_NAME} ${GDAL_LIBRARY})

106
modules/odm_dem2points/src/CmdLineParser.h
Wyświetl plik

@ -1,106 +0,0 @@
/*
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
All rights reserved.
Redistribution and use in source and binary forms, with or without modification,
are permitted provided that the following conditions are met:
Redistributions of source code must retain the above copyright notice, this list of
conditions and the following disclaimer. Redistributions in binary form must reproduce
the above copyright notice, this list of conditions and the following disclaimer
in the documentation and/or other materials provided with the distribution.
Neither the name of the Johns Hopkins University nor the names of its contributors
may be used to endorse or promote products derived from this software without specific
prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
DAMAGE.
*/
#ifndef CMD_LINE_PARSER_INCLUDED
#define CMD_LINE_PARSER_INCLUDED
#include <stdarg.h>
#include <cstring>
#include <cstdlib>
#include <string>
#include <vector>
#ifdef WIN32
int strcasecmp( const char* c1 , const char* c2 );
#endif // WIN32
class cmdLineReadable
{
public:
bool set;
char *name;
cmdLineReadable( const char *name );
virtual ~cmdLineReadable( void );
virtual int read( char** argv , int argc );
virtual void writeValue( char* str ) const;
};
template< class Type > void cmdLineWriteValue( Type t , char* str );
template< class Type > void cmdLineCleanUp( Type* t );
template< class Type > Type cmdLineInitialize( void );
template< class Type > Type cmdLineCopy( Type t );
template< class Type > Type cmdLineStringToType( const char* str );
template< class Type >
class cmdLineParameter : public cmdLineReadable
{
public:
Type value;
cmdLineParameter( const char *name );
cmdLineParameter( const char *name , Type v );
~cmdLineParameter( void );
int read( char** argv , int argc );
void writeValue( char* str ) const;
bool expectsArg( void ) const { return true; }
};
template< class Type , int Dim >
class cmdLineParameterArray : public cmdLineReadable
{
public:
Type values[Dim];
cmdLineParameterArray( const char *name, const Type* v=NULL );
~cmdLineParameterArray( void );
int read( char** argv , int argc );
void writeValue( char* str ) const;
bool expectsArg( void ) const { return true; }
};
template< class Type >
class cmdLineParameters : public cmdLineReadable
{
public:
int count;
Type *values;
cmdLineParameters( const char* name );
~cmdLineParameters( void );
int read( char** argv , int argc );
void writeValue( char* str ) const;
bool expectsArg( void ) const { return true; }
};
void cmdLineParse( int argc , char **argv, cmdLineReadable** params );
char* FileExtension( char* fileName );
char* LocalFileName( char* fileName );
char* DirectoryName( char* fileName );
char* GetFileExtension( const char* fileName );
char* GetLocalFileName( const char* fileName );
char** ReadWords( const char* fileName , int& cnt );
#include "CmdLineParser.inl"
#endif // CMD_LINE_PARSER_INCLUDED

300
modules/odm_dem2points/src/CmdLineParser.inl
Wyświetl plik

@ -1,300 +0,0 @@
/* -*- C++ -*-
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
All rights reserved.
Redistribution and use in source and binary forms, with or without modification,
are permitted provided that the following conditions are met:
Redistributions of source code must retain the above copyright notice, this list of
conditions and the following disclaimer. Redistributions in binary form must reproduce
the above copyright notice, this list of conditions and the following disclaimer
in the documentation and/or other materials provided with the distribution.
Neither the name of the Johns Hopkins University nor the names of its contributors
may be used to endorse or promote products derived from this software without specific
prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
DAMAGE.
*/
#include <cassert>
#include <string.h>
#if defined( WIN32 ) || defined( _WIN64 )
inline int strcasecmp( const char* c1 , const char* c2 ){ return _stricmp( c1 , c2 ); }
#endif // WIN32 || _WIN64
template< > void cmdLineCleanUp< int >( int* t ){ *t = 0; }
template< > void cmdLineCleanUp< float >( float* t ){ *t = 0; }
template< > void cmdLineCleanUp< double >( double* t ){ *t = 0; }
template< > void cmdLineCleanUp< char* >( char** t ){ if( *t ) free( *t ) ; *t = NULL; }
template< > int cmdLineInitialize< int >( void ){ return 0; }
template< > float cmdLineInitialize< float >( void ){ return 0.f; }
template< > double cmdLineInitialize< double >( void ){ return 0.; }
template< > char* cmdLineInitialize< char* >( void ){ return NULL; }
template< > void cmdLineWriteValue< int >( int t , char* str ){ sprintf( str , "%d" , t ); }
template< > void cmdLineWriteValue< float >( float t , char* str ){ sprintf( str , "%f" , t ); }
template< > void cmdLineWriteValue< double >( double t , char* str ){ sprintf( str , "%f" , t ); }
template< > void cmdLineWriteValue< char* >( char* t , char* str ){ if( t ) sprintf( str , "%s" , t ) ; else str[0]=0; }
template< > int cmdLineCopy( int t ){ return t; }
template< > float cmdLineCopy( float t ){ return t; }
template< > double cmdLineCopy( double t ){ return t; }
#if defined( WIN32 ) || defined( _WIN64 )
template< > char* cmdLineCopy( char* t ){ return _strdup( t ); }
#else // !WIN32 && !_WIN64
template< > char* cmdLineCopy( char* t ){ return strdup( t ); }
#endif // WIN32 || _WIN64
template< > int cmdLineStringToType( const char* str ){ return atoi( str ); }
template< > float cmdLineStringToType( const char* str ){ return float( atof( str ) ); }
template< > double cmdLineStringToType( const char* str ){ return double( atof( str ) ); }
#if defined( WIN32 ) || defined( _WIN64 )
template< > char* cmdLineStringToType( const char* str ){ return _strdup( str ); }
#else // !WIN32 && !_WIN64
template< > char* cmdLineStringToType( const char* str ){ return strdup( str ); }
#endif // WIN32 || _WIN64
/////////////////////
// cmdLineReadable //
/////////////////////
#if defined( WIN32 ) || defined( _WIN64 )
inline cmdLineReadable::cmdLineReadable( const char *name ) : set(false) { this->name = _strdup( name ); }
#else // !WIN32 && !_WIN64
inline cmdLineReadable::cmdLineReadable( const char *name ) : set(false) { this->name = strdup( name ); }
#endif // WIN32 || _WIN64
inline cmdLineReadable::~cmdLineReadable( void ){ if( name ) free( name ) ; name = NULL; }
inline int cmdLineReadable::read( char** , int ){ set = true ; return 0; }
inline void cmdLineReadable::writeValue( char* str ) const { str[0] = 0; }
//////////////////////
// cmdLineParameter //
//////////////////////
template< class Type > cmdLineParameter< Type >::~cmdLineParameter( void ) { cmdLineCleanUp( &value ); }
template< class Type > cmdLineParameter< Type >::cmdLineParameter( const char *name ) : cmdLineReadable( name ){ value = cmdLineInitialize< Type >(); }
template< class Type > cmdLineParameter< Type >::cmdLineParameter( const char *name , Type v ) : cmdLineReadable( name ){ value = cmdLineCopy< Type >( v ); }
template< class Type >
int cmdLineParameter< Type >::read( char** argv , int argc )
{
if( argc>0 )
{
cmdLineCleanUp< Type >( &value ) , value = cmdLineStringToType< Type >( argv[0] );
set = true;
return 1;
}
else return 0;
}
template< class Type >
void cmdLineParameter< Type >::writeValue( char* str ) const { cmdLineWriteValue< Type >( value , str ); }
///////////////////////////
// cmdLineParameterArray //
///////////////////////////
template< class Type , int Dim >
cmdLineParameterArray< Type , Dim >::cmdLineParameterArray( const char *name , const Type* v ) : cmdLineReadable( name )
{
if( v ) for( int i=0 ; i<Dim ; i++ ) values[i] = cmdLineCopy< Type >( v[i] );
else for( int i=0 ; i<Dim ; i++ ) values[i] = cmdLineInitialize< Type >();
}
template< class Type , int Dim >
cmdLineParameterArray< Type , Dim >::~cmdLineParameterArray( void ){ for( int i=0 ; i<Dim ; i++ ) cmdLineCleanUp< Type >( values+i ); }
template< class Type , int Dim >
int cmdLineParameterArray< Type , Dim >::read( char** argv , int argc )
{
if( argc>=Dim )
{
for( int i=0 ; i<Dim ; i++ ) cmdLineCleanUp< Type >( values+i ) , values[i] = cmdLineStringToType< Type >( argv[i] );
set = true;
return Dim;
}
else return 0;
}
template< class Type , int Dim >
void cmdLineParameterArray< Type , Dim >::writeValue( char* str ) const
{
char* temp=str;
for( int i=0 ; i<Dim ; i++ )
{
cmdLineWriteValue< Type >( values[i] , temp );
temp = str+strlen( str );
}
}
///////////////////////
// cmdLineParameters //
///////////////////////
template< class Type >
cmdLineParameters< Type >::cmdLineParameters( const char* name ) : cmdLineReadable( name ) , values(NULL) , count(0) { }
template< class Type >
cmdLineParameters< Type >::~cmdLineParameters( void )
{
if( values ) delete[] values;
values = NULL;
count = 0;
}
template< class Type >
int cmdLineParameters< Type >::read( char** argv , int argc )
{
if( values ) delete[] values;
values = NULL;
if( argc>0 )
{
count = atoi(argv[0]);
if( count <= 0 || argc <= count ) return 1;
values = new Type[count];
if( !values ) return 0;
for( int i=0 ; i<count ; i++ ) values[i] = cmdLineStringToType< Type >( argv[i+1] );
set = true;
return count+1;
}
else return 0;
}
template< class Type >
void cmdLineParameters< Type >::writeValue( char* str ) const
{
char* temp=str;
for( int i=0 ; i<count ; i++ )
{
cmdLineWriteValue< Type >( values[i] , temp );
temp = str+strlen( str );
}
}
inline char* FileExtension( char* fileName )
{
char* temp = fileName;
for( unsigned int i=0 ; i<strlen(fileName) ; i++ ) if( fileName[i]=='.' ) temp = &fileName[i+1];
return temp;
}
inline char* GetFileExtension( const char* fileName )
{
char* fileNameCopy;
char* ext=NULL;
char* temp;
fileNameCopy=new char[strlen(fileName)+1];
assert(fileNameCopy);
strcpy(fileNameCopy,fileName);
temp=strtok(fileNameCopy,".");
while(temp!=NULL)
{
if(ext!=NULL){delete[] ext;}
ext=new char[strlen(temp)+1];
assert(ext);
strcpy(ext,temp);
temp=strtok(NULL,".");
}
delete[] fileNameCopy;
return ext;
}
inline char* GetLocalFileName( const char* fileName )
{
char* fileNameCopy;
char* name=NULL;
char* temp;
fileNameCopy=new char[strlen(fileName)+1];
assert(fileNameCopy);
strcpy(fileNameCopy,fileName);
temp=strtok(fileNameCopy,"\\");
while(temp!=NULL){
if(name!=NULL){delete[] name;}
name=new char[strlen(temp)+1];
assert(name);
strcpy(name,temp);
temp=strtok(NULL,"\\");
}
delete[] fileNameCopy;
return name;
}
inline char* LocalFileName( char* fileName )
{
char* temp = fileName;
for( int i=0 ; i<(int)strlen(fileName) ; i++ ) if( fileName[i] =='\\' ) temp = &fileName[i+1];
return temp;
}
inline char* DirectoryName( char* fileName )
{
for( int i=int( strlen(fileName) )-1 ; i>=0 ; i-- )
if( fileName[i] =='\\' )
{
fileName[i] = 0;
return fileName;
}
fileName[0] = 0;
return fileName;
}
inline void cmdLineParse( int argc , char **argv , cmdLineReadable** params )
{
while( argc>0 )
{
if( argv[0][0]=='-' )
{
cmdLineReadable* readable=NULL;
for( int i=0 ; params[i]!=NULL && readable==NULL ; i++ ) if( !strcasecmp( params[i]->name , argv[0]+1 ) ) readable = params[i];
if( readable )
{
int j = readable->read( argv+1 , argc-1 );
argv += j , argc -= j;
}
else
{
fprintf( stderr , "[WARNING] Invalid option: %s\n" , argv[0] );
for( int i=0 ; params[i]!=NULL ; i++ ) printf( "\t-%s\n" , params[i]->name );
}
}
else fprintf( stderr , "[WARNING] Parameter name should be of the form -<name>: %s\n" , argv[0] );
++argv , --argc;
}
}
inline char** ReadWords(const char* fileName,int& cnt)
{
char** names;
char temp[500];
FILE* fp;
fp=fopen(fileName,"r");
if(!fp){return NULL;}
cnt=0;
while(fscanf(fp," %s ",temp)==1){cnt++;}
fclose(fp);
names=new char*[cnt];
if(!names){return NULL;}
fp=fopen(fileName,"r");
if(!fp){
delete[] names;
cnt=0;
return NULL;
}
cnt=0;
while(fscanf(fp," %s ",temp)==1){
names[cnt]=new char[strlen(temp)+1];
if(!names){
for(int j=0;j<cnt;j++){delete[] names[j];}
delete[] names;
cnt=0;
fclose(fp);
return NULL;
}
strcpy(names[cnt],temp);
cnt++;
}
fclose(fp);
return names;
}

33
modules/odm_dem2points/src/Logger.h
Wyświetl plik

@ -1,33 +0,0 @@
#include <cstdio>
#include <cstdarg>
#include "CmdLineParser.h"
struct Logger{
bool verbose;
const char* outputFile;
Logger(){
this->verbose = false;
this->outputFile = NULL;
}
void operator() ( const char* format , ... )
{
if( outputFile )
{
FILE* fp = fopen( outputFile , "a" );
va_list args;
va_start( args , format );
vfprintf( fp , format , args );
fclose( fp );
va_end( args );
}
if( verbose )
{
va_list args;
va_start( args , format );
vprintf( format , args );
va_end( args );
}
}
};

234
modules/odm_dem2points/src/main.cpp
Wyświetl plik

@ -1,234 +0,0 @@
#include <iostream>
#include <string>
#include <fstream>
#include "CmdLineParser.h"
#include "Logger.h"
#include "gdal_priv.h"
#include "cpl_conv.h" // for CPLMalloc()
Logger logWriter;
typedef struct Point2D{
double x;
double y;
Point2D(): x(0.0), y(0.0){}
Point2D(double x, double y): x(x), y(y){}
} Point2D;
typedef struct BoundingBox{
Point2D max;
Point2D min;
BoundingBox(): max(Point2D()), min(Point2D()){}
BoundingBox(Point2D min, Point2D max): max(max), min(min){}
} BoundingBox;
struct PlyPoint{
float x;
float y;
float z;
float nx;
float ny;
float nz;
} p;
size_t psize = sizeof(float) * 6;
float *rasterData;
int arr_width, arr_height;
#define IS_BIG_ENDIAN (*(uint16_t *)"\0\xff" < 0x100)
cmdLineParameter< char* >
InputFile( "inputFile" ) ,
OutputFile( "outputFile" );
cmdLineParameter< float >
SkirtHeightThreshold( "skirtHeightThreshold" ),
SkirtIncrements("skirtIncrements"),
SkirtHeightCap( "skirtHeightCap");
cmdLineReadable
Verbose( "verbose" );
cmdLineReadable* params[] = {
&InputFile , &OutputFile , &SkirtHeightThreshold, &SkirtIncrements, &SkirtHeightCap, &Verbose ,
NULL
};
void help(char *ex){
std::cout << "Usage: " << ex << std::endl
<< "\t -" << InputFile.name << " <input DSM raster>" << std::endl
<< "\t -" << OutputFile.name << " <output PLY points>" << std::endl
<< "\t [-" << SkirtHeightThreshold.name << " <Height threshold between cells that triggers the creation of a skirt>]" << std::endl
<< "\t [-" << SkirtIncrements.name << " <Skirt height increments when adding a new skirt>]" << std::endl
<< "\t [-" << SkirtHeightCap.name << " <Height cap that blocks the creation of a skirt>]" << std::endl
<< "\t [-" << Verbose.name << "]" << std::endl;
exit(EXIT_FAILURE);
}
void logArgs(cmdLineReadable* params[], Logger& logWriter){
logWriter("Running with parameters:\n");
char str[1024];
for( int i=0 ; params[i] ; i++ ){
if( params[i]->set ){
params[i]->writeValue( str );
if( strlen( str ) ) logWriter( "\t--%s %s\n" , params[i]->name , str );
else logWriter( "\t--%s\n" , params[i]->name );
}
}
}
Point2D geoLoc(float xloc, float yloc, double *affine){
return Point2D(affine[0] + xloc*affine[1] + yloc*affine[2], affine[3] + xloc*affine[4] + yloc*affine[5]);
}
BoundingBox getExtent(GDALDataset *dataset){
double affine[6];
dataset->GetGeoTransform(affine);
return BoundingBox(geoLoc(0, dataset->GetRasterYSize(), affine), geoLoc(dataset->GetRasterXSize(), 0, affine));
}
int countSkirts(int nx, int ny){
float neighbor_z = rasterData[ny * arr_width + nx];
float current_z = p.z;
int result = 0;
if (current_z - neighbor_z > SkirtHeightThreshold.value && current_z - neighbor_z < SkirtHeightCap.value){
while (current_z > neighbor_z){
current_z -= SkirtIncrements.value;
result++;
}
result++;
}
return result;
}
void writeSkirts(std::ofstream &f, int nx, int ny){
float neighbor_z = rasterData[ny * arr_width + nx];
float z_val = p.z;
if (p.z - neighbor_z > SkirtHeightThreshold.value && p.z - neighbor_z < SkirtHeightCap.value){
while (p.z > neighbor_z){
p.z -= SkirtIncrements.value;
f.write(reinterpret_cast<char*>(&p), psize);
}
p.z = neighbor_z;
f.write(reinterpret_cast<char*>(&p), psize);
}
// Restore original p.z value
p.z = z_val;
}
int main(int argc, char **argv) {
cmdLineParse( argc-1 , &argv[1] , params );
if( !InputFile.set || !OutputFile.set ) help(argv[0]);
if (!SkirtHeightThreshold.set) SkirtHeightThreshold.value = 1.5f;
if (!SkirtIncrements.set) SkirtIncrements.value = 0.1f;
if (!SkirtHeightCap.set) SkirtHeightCap.value = 100.0f;
logWriter.verbose = Verbose.set;
// logWriter.outputFile = "odm_dem2points_log.txt";
logArgs(params, logWriter);
GDALDataset *dataset;
GDALAllRegister();
dataset = (GDALDataset *) GDALOpen( InputFile.value, GA_ReadOnly );
if( dataset != NULL )
{
arr_width = dataset->GetRasterXSize();
arr_height = dataset->GetRasterYSize();
logWriter("Raster Size is %dx%d\n", arr_width, arr_height);
BoundingBox extent = getExtent(dataset);
logWriter("Extent is (%f, %f), (%f, %f)\n", extent.min.x, extent.max.x, extent.min.y, extent.max.y);
float ext_width = extent.max.x - extent.min.x;
float ext_height = extent.max.y - extent.min.y;
int vertex_count = (arr_height - 4) * (arr_width - 4);
int skirt_points = 0;
GDALRasterBand *band = dataset->GetRasterBand(1);
rasterData = new float[arr_width * arr_height];
if (band->RasterIO( GF_Read, 0, 0, arr_width, arr_height,
rasterData, arr_width, arr_height, GDT_Float32, 0, 0 ) == CE_Failure){
std::cerr << "Cannot access raster data" << std::endl;
exit(1);
}
// On the first pass we calculate the number of skirts
// on the second pass we sample and write the points
logWriter("Calculating skirts... ");
for (int y = 2; y < arr_height - 2; y++){
for (int x = 2; x < arr_width - 2; x++){
p.z = rasterData[y * arr_width + x];
skirt_points += countSkirts(x, y + 1);
skirt_points += countSkirts(x, y - 1);
skirt_points += countSkirts(x + 1, y);
skirt_points += countSkirts(x - 1, y);
}
}
logWriter("%d vertices will be added\n", skirt_points);
logWriter("Total vertices: %d\n", (skirt_points + vertex_count));
logWriter("Sampling and writing to file...");
// Starting writing ply file
std::ofstream f (OutputFile.value);
f << "ply" << std::endl;
if (IS_BIG_ENDIAN){
f << "format binary_big_endian 1.0" << std::endl;
}else{
f << "format binary_little_endian 1.0" << std::endl;
}
f << "element vertex " << (vertex_count + skirt_points) << std::endl
<< "property float x" << std::endl
<< "property float y" << std::endl
<< "property float z" << std::endl
<< "property float nx" << std::endl
<< "property float ny" << std::endl
<< "property float nz" << std::endl
<< "end_header" << std::endl;
// Normals are always pointing up
p.nx = 0;
p.ny = 0;
p.nz = 1;
for (int y = 2; y < arr_height - 2; y++){
for (int x = 2; x < arr_width - 2; x++){
p.z = rasterData[y * arr_width + x];
p.x = extent.min.x + (static_cast<float>(x) / static_cast<float>(arr_width)) * ext_width;
p.y = extent.max.y - (static_cast<float>(y) / static_cast<float>(arr_height)) * ext_height;
f.write(reinterpret_cast<char*>(&p), psize);
writeSkirts(f, x, y + 1);
writeSkirts(f, x, y - 1);
writeSkirts(f, x + 1, y);
writeSkirts(f, x - 1, y);
}
}
logWriter(" done!\n");
f.close();
GDALClose(dataset);
}else{
std::cerr << "Cannot open " << InputFile.value << std::endl;
}
}

3
opendm/context.py
Wyświetl plik

@ -28,6 +28,9 @@ makescene_path = os.path.join(superbuild_path, 'src', 'elibs', 'mve', 'apps', 'm
smvs_path = os.path.join(superbuild_path, 'src', 'elibs', 'smvs', 'app', 'smvsrecon')
poisson_recon_path = os.path.join(superbuild_path, 'src', 'PoissonRecon', 'Bin', 'Linux', 'PoissonRecon')
dem2mesh_path = os.path.join(superbuild_path, 'src', 'dem2mesh', 'dem2mesh')
dem2points_path = os.path.join(superbuild_path, 'src', 'dem2points', 'dem2points')
# define mvstex path
mvstex_path = os.path.join(superbuild_path, "install/bin/texrecon")

11
opendm/dem/commands.py
Wyświetl plik

@ -127,10 +127,19 @@ def gap_fill(filenames, fout):
return_indices=True)
arr = arr[tuple(indices)]
# Median filter
# Median filter (careful, changing the value 5 might require tweaking)
# the lines below. There's another numpy function that takes care of
# these edge cases, but it's slower.
from scipy import signal
arr = signal.medfilt(arr, 5)
# Fill corner points with nearest value
if arr.shape >= (4, 4):
arr[0][:2] = arr[1][0] = arr[1][1]
arr[0][-2:] = arr[1][-1] = arr[2][-1]
arr[-1][:2] = arr[-2][0] = arr[-2][1]
arr[-1][-2:] = arr[-2][-1] = arr[-2][-2]
# write output
imgout = gippy.GeoImage.create_from(imgs[0], fout)
imgout.set_nodata(nodata)

18
opendm/gsd.py
Wyświetl plik

@ -4,6 +4,22 @@ import numpy as np
from repoze.lru import lru_cache
from opendm import log
def rounded_gsd(reconstruction_json, default_value=None, ndigits=0, ignore_gsd=False):
"""
:param reconstruction_json path to OpenSfM's reconstruction.json
:return GSD value rounded. If GSD cannot be computed, or ignore_gsd is set, it returns a default value.
"""
if ignore_gsd:
return default_value
gsd = opensfm_reconstruction_average_gsd(reconstruction_json)
if gsd is not None:
return round(gsd, ndigits)
else:
return default_value
def image_scale_factor(target_resolution, reconstruction_json, gsd_error_estimate = 0.1):
"""
:param target_resolution resolution the user wants have in cm / pixel
@ -56,7 +72,7 @@ def opensfm_reconstruction_average_gsd(reconstruction_json):
a GSD estimate cannot be compute
"""
if not os.path.isfile(reconstruction_json):
raise FileNotFoundError(reconstruction_json + " does not exist.")
raise IOError(reconstruction_json + " does not exist.")
with open(reconstruction_json) as f:
data = json.load(f)

91
opendm/mesh.py
Wyświetl plik

@ -8,7 +8,7 @@ from opendm import context
from scipy import signal, ndimage
import numpy as np
def create_25dmesh(inPointCloud, outMesh, dsm_resolution=0.05, depth=8, samples=1, maxVertexCount=100000, verbose=False, max_workers=None):
def create_25dmesh(inPointCloud, outMesh, dsm_radius=0.07, dsm_resolution=0.05, depth=8, samples=1, maxVertexCount=100000, verbose=False, max_workers=None, method='gridded'):
# Create DSM from point cloud
# Create temporary directory
@ -19,7 +19,7 @@ def create_25dmesh(inPointCloud, outMesh, dsm_resolution=0.05, depth=8, samples=
os.mkdir(tmp_directory)
log.ODM_INFO('Created temporary directory: %s' % tmp_directory)
radius_steps = [dsm_resolution * math.sqrt(2)]
radius_steps = [dsm_radius]
log.ODM_INFO('Creating DSM for 2.5D mesh')
@ -30,17 +30,22 @@ def create_25dmesh(inPointCloud, outMesh, dsm_resolution=0.05, depth=8, samples=
gapfill=True,
outdir=tmp_directory,
resolution=dsm_resolution,
products=['idw'],
products=['max'],
verbose=verbose,
max_workers=max_workers
)
dsm_points = dem_to_points(os.path.join(tmp_directory, 'mesh_dsm.tif'), os.path.join(tmp_directory, 'dsm_points.ply'), verbose)
mesh = screened_poisson_reconstruction(dsm_points, outMesh, depth=depth,
if method == 'gridded':
mesh = dem_to_mesh_gridded(os.path.join(tmp_directory, 'mesh_dsm.tif'), outMesh, maxVertexCount, verbose)
elif method == 'poisson':
dsm_points = dem_to_points(os.path.join(tmp_directory, 'mesh_dsm.tif'), os.path.join(tmp_directory, 'dsm_points.ply'), verbose)
mesh = screened_poisson_reconstruction(dsm_points, outMesh, depth=depth,
samples=samples,
maxVertexCount=maxVertexCount,
threads=max_workers,
verbose=verbose)
else:
raise 'Not a valid method: ' + method
# Cleanup tmp
if os.path.exists(tmp_directory):
@ -48,17 +53,18 @@ def create_25dmesh(inPointCloud, outMesh, dsm_resolution=0.05, depth=8, samples=
return mesh
def dem_to_points(inGeotiff, outPointCloud, verbose=False):
log.ODM_INFO('Sampling points from DSM: %s' % inGeotiff)
kwargs = {
'bin': context.odm_modules_path,
'bin': context.dem2points_path,
'outfile': outPointCloud,
'infile': inGeotiff,
'verbose': '-verbose' if verbose else ''
}
system.run('{bin}/odm_dem2points -inputFile {infile} '
system.run('{bin} -inputFile {infile} '
'-outputFile {outfile} '
'-skirtHeightThreshold 1.5 '
'-skirtIncrements 0.2 '
@ -67,67 +73,24 @@ def dem_to_points(inGeotiff, outPointCloud, verbose=False):
return outPointCloud
# Old Python implementation of dem_to_points
# def dem_to_points(inGeotiff, outPointCloud):
# log.ODM_INFO('Sampling points from DSM: %s' % inGeotiff)
# image = GeoImage.open([inGeotiff], bandnames=['z'], nodata=-9999)
# arr = image['z'].read_raw()
def dem_to_mesh_gridded(inGeotiff, outPointCloud, maxVertexCount, verbose=False):
log.ODM_INFO('Creating mesh from DSM: %s' % inGeotiff)
# mem_file = BytesIO()
kwargs = {
'bin': context.dem2mesh_path,
'outfile': outPointCloud,
'infile': inGeotiff,
'maxVertexCount': maxVertexCount,
'verbose': '-verbose' if verbose else ''
}
# xmin, xmax, ymin, ymax = image.extent().x0(), image.extent().x1(), image.extent().y0(), image.extent().y1()
# ext_width, ext_height = xmax - xmin, ymax - ymin
# arr_height, arr_width = arr.shape
# vertex_count = (arr_height - 4) * (arr_width - 4)
# skirt_points = 0
system.run('{bin} -inputFile {infile} '
'-outputFile {outfile} '
'-maxVertexCount {maxVertexCount} '
' {verbose} '.format(**kwargs))
# skirt_height_threshold = 1 # meter
# skirt_increments = 0.1
# for x in range(2, arr_width - 2):
# for y in range(2, arr_height - 2):
# z = arr[y][x]
# tx = xmin + (float(x) / float(arr_width)) * ext_width
# ty = ymax - (float(y) / float(arr_height)) * ext_height
# mem_file.write(struct.pack('ffffff', tx, ty, z, 0, 0, 1))
# # Skirting
# for (nx, ny) in ((x, y + 1), (x, y - 1), (x + 1, y), (x - 1, y)):
# current_z = z
# neighbor_z = arr[ny][nx]
# if current_z - neighbor_z > skirt_height_threshold:
# while current_z > neighbor_z:
# current_z -= skirt_increments
# mem_file.write(struct.pack('ffffff', tx, ty, current_z, 0, 0, 1))
# skirt_points += 1
# mem_file.write(struct.pack('ffffff', tx, ty, neighbor_z, 0, 0, 1))
# skirt_points += 1
# log.ODM_INFO("Points count: %s (%s samples, %s skirts)", vertex_count + skirt_points, vertex_count, skirt_points)
# log.ODM_INFO('Writing points...')
# mem_file.seek(0)
# with open(outPointCloud, "wb") as f:
# f.write("ply\n")
# f.write("format binary_%s_endian 1.0\n" % sys.byteorder)
# f.write("element vertex %s\n" % (vertex_count + skirt_points))
# f.write("property float x\n")
# f.write("property float y\n")
# f.write("property float z\n")
# f.write("property float nx\n")
# f.write("property float ny\n")
# f.write("property float nz\n")
# f.write("end_header\n")
# shutil.copyfileobj(mem_file, f)
# mem_file.close()
# log.ODM_INFO('Wrote points to: %s' % outPointCloud)
# return outPointCloud
return outPointCloud
def screened_poisson_reconstruction(inPointCloud, outMesh, depth = 8, samples = 1, maxVertexCount=100000, pointWeight=4, threads=context.num_cores, verbose=False):

28
scripts/odm_meshing.py
Wyświetl plik

@ -1,4 +1,4 @@
import ecto, os
import ecto, os, math
from opendm import log
from opendm import io
@ -79,23 +79,33 @@ class ODMeshingCell(ecto.Cell):
if not io.file_exists(tree.odm_25dmesh) or rerun_cell:
log.ODM_DEBUG('Writing ODM 2.5D Mesh file in: %s' % tree.odm_25dmesh)
dsm_resolution = gsd.cap_resolution(args.orthophoto_resolution, tree.opensfm_reconstruction, ignore_gsd=args.ignore_gsd) / 100.0
# Create reference DSM at half ortho resolution
dsm_resolution *= 2
ortho_resolution = gsd.cap_resolution(args.orthophoto_resolution, tree.opensfm_reconstruction, ignore_gsd=args.ignore_gsd) / 100.0
dsm_multiplier = max(1.0, gsd.rounded_gsd(tree.opensfm_reconstruction, default_value=4, ndigits=3, ignore_gsd=args.ignore_gsd))
# A good DSM size depends on the flight altitude.
# Flights at low altitude need more details (higher resolution)
# Flights at higher altitude benefit from smoother surfaces (lower resolution)
dsm_resolution = ortho_resolution * dsm_multiplier
dsm_radius = dsm_resolution * math.sqrt(2)
# Sparse point clouds benefits from using
# a larger resolution value (more radius interolation, less holes)
# a larger radius interolation --> less holes
if args.fast_orthophoto:
dsm_resolution *= 2
dsm_radius *= 2
log.ODM_DEBUG('ODM 2.5D DSM resolution: %s' % dsm_resolution)
mesh.create_25dmesh(infile, tree.odm_25dmesh,
dsm_resolution=dsm_resolution,
dsm_radius=dsm_radius,
dsm_resolution=dsm_resolution,
depth=self.params.oct_tree,
maxVertexCount=self.params.max_vertex,
samples=self.params.samples,
verbose=self.params.verbose,
max_workers=args.max_concurrency)
max_workers=args.max_concurrency,
method='poisson' if args.fast_orthophoto else 'gridded')
else:
log.ODM_WARNING('Found a valid ODM 2.5D Mesh file in: %s' %
tree.odm_25dmesh)

5
scripts/smvs.py
Wyświetl plik

@ -48,10 +48,13 @@ class ODMSmvsCell(ecto.Cell):
# check if reconstruction was done before
if not io.file_exists(tree.smvs_model) or rerun_cell:
# cleanup if a rerun
if io.dir_exists(tree.mve_path) and rerun_cell:
shutil.rmtree(tree.mve_path)
# make bundle directory
if not io.file_exists(tree.mve_bundle):
system.mkdir_p(tree.mve_path) #TODO: check permissions/what happens when rerun
system.mkdir_p(tree.mve_path)
system.mkdir_p(io.join_paths(tree.mve_path, 'bundle'))
io.copy(tree.opensfm_image_list, tree.mve_image_list)
io.copy(tree.opensfm_bundle, tree.mve_bundle)