1
0
Fork 0

Add a tool to chop CGIAR (http://srtm.csi.cgiar.org/) srtm improved files

This commit is contained in:
Frederic Bouvier 2009-02-23 09:14:59 +01:00 committed by Ralf Gerlich
parent 1bca4c11e2
commit 00b8658b99
11 changed files with 827 additions and 153 deletions

View file

@ -104,6 +104,9 @@
<File
RelativePath="..\..\..\src\Lib\Hgt\hgt.cxx">
</File>
<File
RelativePath="..\..\..\src\Lib\Hgt\srtmbase.cxx">
</File>
</Filter>
<Filter
Name="Header Files"
@ -112,6 +115,9 @@
<File
RelativePath="..\..\..\src\Lib\Hgt\hgt.hxx">
</File>
<File
RelativePath="..\..\..\src\Lib\Hgt\srtmbase.hxx">
</File>
</Filter>
<Filter
Name="Resource Files"

2
projects/VC7.1/srtmchop/.gitignore vendored Normal file
View file

@ -0,0 +1,2 @@
Debug
Release

View file

@ -0,0 +1,140 @@
<?xml version="1.0" encoding="Windows-1252"?>
<VisualStudioProject
ProjectType="Visual C++"
Version="7.10"
Name="srtmchop"
ProjectGUID="{7490AAD8-D2A4-4F0D-8150-68DCE07F3116}"
Keyword="Win32Proj">
<Platforms>
<Platform
Name="Win32"/>
</Platforms>
<Configurations>
<Configuration
Name="Debug|Win32"
OutputDirectory="..\bin\Debug"
IntermediateDirectory="Debug"
ConfigurationType="1"
CharacterSet="2">
<Tool
Name="VCCLCompilerTool"
Optimization="0"
AdditionalIncludeDirectories="..\..\..\..\SimGear.cs;..\..\..\..;..\..\..\..\3rdParty\include;..\..\..\..\install\msvc71\OpenSceneGraph\include;..\..\..\src\Lib;..\..\..\src\BuildTiles;..\..\..\..\gpc232;..\..\..\src"
PreprocessorDefinitions="WIN32;_DEBUG;_CONSOLE;NOMINMAX"
MinimalRebuild="TRUE"
BasicRuntimeChecks="3"
RuntimeLibrary="3"
RuntimeTypeInfo="TRUE"
UsePrecompiledHeader="0"
WarningLevel="3"
Detect64BitPortabilityProblems="TRUE"
DebugInformationFormat="4"/>
<Tool
Name="VCCustomBuildTool"/>
<Tool
Name="VCLinkerTool"
AdditionalDependencies="ul_d.lib sg_d.lib zlibd.lib libtiffd.lib"
OutputFile="$(OutDir)/srtmchop.exe"
LinkIncremental="2"
AdditionalLibraryDirectories="..\..\..\..\plib;..\..\..\..\3rdParty\lib"
GenerateDebugInformation="TRUE"
ProgramDatabaseFile="$(OutDir)/srtmchop.pdb"
SubSystem="1"
TargetMachine="1"/>
<Tool
Name="VCMIDLTool"/>
<Tool
Name="VCPostBuildEventTool"/>
<Tool
Name="VCPreBuildEventTool"/>
<Tool
Name="VCPreLinkEventTool"/>
<Tool
Name="VCResourceCompilerTool"/>
<Tool
Name="VCWebServiceProxyGeneratorTool"/>
<Tool
Name="VCXMLDataGeneratorTool"/>
<Tool
Name="VCWebDeploymentTool"/>
<Tool
Name="VCManagedWrapperGeneratorTool"/>
<Tool
Name="VCAuxiliaryManagedWrapperGeneratorTool"/>
</Configuration>
<Configuration
Name="Release|Win32"
OutputDirectory="..\bin\Release"
IntermediateDirectory="Release"
ConfigurationType="1"
CharacterSet="2">
<Tool
Name="VCCLCompilerTool"
AdditionalIncludeDirectories="..\..\..\..\SimGear.cs;..\..\..\..;..\..\..\..\3rdParty\include;..\..\..\..\install\msvc71\OpenSceneGraph\include;..\..\..\src\Lib;..\..\..\src\BuildTiles;..\..\..\..\gpc232;..\..\..\src"
PreprocessorDefinitions="WIN32;NDEBUG;_CONSOLE;NOMINMAX"
RuntimeLibrary="2"
RuntimeTypeInfo="TRUE"
UsePrecompiledHeader="0"
WarningLevel="3"
Detect64BitPortabilityProblems="TRUE"
DebugInformationFormat="3"/>
<Tool
Name="VCCustomBuildTool"/>
<Tool
Name="VCLinkerTool"
AdditionalDependencies="ul.lib sg.lib zlib.lib libtiff.lib"
OutputFile="$(OutDir)/srtmchop.exe"
LinkIncremental="1"
AdditionalLibraryDirectories="..\..\..\..\plib;..\..\..\..\3rdParty\lib"
GenerateDebugInformation="TRUE"
SubSystem="1"
OptimizeReferences="2"
EnableCOMDATFolding="2"
TargetMachine="1"/>
<Tool
Name="VCMIDLTool"/>
<Tool
Name="VCPostBuildEventTool"/>
<Tool
Name="VCPreBuildEventTool"/>
<Tool
Name="VCPreLinkEventTool"/>
<Tool
Name="VCResourceCompilerTool"/>
<Tool
Name="VCWebServiceProxyGeneratorTool"/>
<Tool
Name="VCXMLDataGeneratorTool"/>
<Tool
Name="VCWebDeploymentTool"/>
<Tool
Name="VCManagedWrapperGeneratorTool"/>
<Tool
Name="VCAuxiliaryManagedWrapperGeneratorTool"/>
</Configuration>
</Configurations>
<References>
</References>
<Files>
<Filter
Name="Source Files"
Filter="cpp;c;cxx;def;odl;idl;hpj;bat;asm;asmx"
UniqueIdentifier="{4FC737F1-C7A5-4376-A066-2A32D752A2FF}">
<File
RelativePath="..\..\..\src\Prep\DemChop\srtmchop.cxx">
</File>
</Filter>
<Filter
Name="Header Files"
Filter="h;hpp;hxx;hm;inl;inc;xsd"
UniqueIdentifier="{93995380-89BD-4b04-88EB-625FBE52EBFB}">
</Filter>
<Filter
Name="Resource Files"
Filter="rc;ico;cur;bmp;dlg;rc2;rct;bin;rgs;gif;jpg;jpeg;jpe;resx"
UniqueIdentifier="{67DA6AB6-F800-4c08-8B7A-83BB121AAD01}">
</Filter>
</Files>
<Globals>
</Globals>
</VisualStudioProject>

View file

@ -233,6 +233,12 @@ Project("{8BC9CEB8-8B4A-11D0-8D11-00A0C91BC942}") = "poly2ogr", "poly2ogr\poly2o
{22540CD3-D3CA-4C86-A773-80AEEE3ACDED} = {22540CD3-D3CA-4C86-A773-80AEEE3ACDED}
EndProjectSection
EndProject
Project("{8BC9CEB8-8B4A-11D0-8D11-00A0C91BC942}") = "srtmchop", "srtmchop\srtmchop.vcproj", "{7490AAD8-D2A4-4F0D-8150-68DCE07F3116}"
ProjectSection(ProjectDependencies) = postProject
{322E9356-2E79-49A6-A7B6-643B35DEB67D} = {322E9356-2E79-49A6-A7B6-643B35DEB67D}
{22540CD3-D3CA-4C86-A773-80AEEE3ACDED} = {22540CD3-D3CA-4C86-A773-80AEEE3ACDED}
EndProjectSection
EndProject
Global
GlobalSection(SolutionConfiguration) = preSolution
Debug = Debug
@ -403,6 +409,10 @@ Global
{74F5BBD4-D669-4348-8B46-11700135CCF1}.Debug.Build.0 = Debug|Win32
{74F5BBD4-D669-4348-8B46-11700135CCF1}.Release.ActiveCfg = Release|Win32
{74F5BBD4-D669-4348-8B46-11700135CCF1}.Release.Build.0 = Release|Win32
{7490AAD8-D2A4-4F0D-8150-68DCE07F3116}.Debug.ActiveCfg = Debug|Win32
{7490AAD8-D2A4-4F0D-8150-68DCE07F3116}.Debug.Build.0 = Debug|Win32
{7490AAD8-D2A4-4F0D-8150-68DCE07F3116}.Release.ActiveCfg = Release|Win32
{7490AAD8-D2A4-4F0D-8150-68DCE07F3116}.Release.Build.0 = Release|Win32
EndGlobalSection
GlobalSection(ExtensibilityGlobals) = postSolution
EndGlobalSection

View file

@ -1,6 +1,6 @@
noinst_LIBRARIES = libHGT.a
libHGT_a_SOURCES = hgt.cxx hgt.hxx
libHGT_a_SOURCES = hgt.cxx hgt.hxx srtmbase.cxx srtmbase.hxx
noinst_PROGRAMS = testhgt

View file

@ -183,124 +183,9 @@ TGHgt::load( ) {
}
// write out the area of data covered by the specified bucket. Data
// is written out column by column starting at the lower left hand
// corner.
bool
TGHgt::write_area( const string& root, SGBucket& b ) {
// calculate some boundaries
double min_x = ( b.get_center_lon() - 0.5 * b.get_width() ) * 3600.0;
double max_x = ( b.get_center_lon() + 0.5 * b.get_width() ) * 3600.0;
double min_y = ( b.get_center_lat() - 0.5 * b.get_height() ) * 3600.0;
double max_y = ( b.get_center_lat() + 0.5 * b.get_height() ) * 3600.0;
cout << b << endl;
cout << "width = " << b.get_width() << " height = " << b.get_height()
<< endl;
cout << "min = " << min_x << "," << min_y
<< " max = " << max_x << "," << max_y << endl;
int start_x = (int)((min_x - originx) / col_step);
int span_x = (int)(b.get_width() * 3600.0 / col_step);
int start_y = (int)((min_y - originy) / row_step);
int span_y = (int)(b.get_height() * 3600.0 / row_step);
cout << "start_x = " << start_x << " span_x = " << span_x << endl;
cout << "start_y = " << start_y << " span_y = " << span_y << endl;
// Do a simple sanity checking. But, please, please be nice to
// this write_area() routine and feed it buckets that coincide
// well with the underlying grid structure and spacing.
if ( ( min_x < originx )
|| ( max_x > originx + cols * col_step )
|| ( min_y < originy )
|| ( max_y > originy + rows * row_step ) ) {
cout << " ERROR: bucket at least partially outside HGT data range!" <<
endl;
return false;
}
// If the area is all ocean, skip it.
if ( !has_non_zero_elev(start_x, span_x, start_y, span_y) ) {
cout << "Tile is all zero elevation: skipping" << endl;
return false;
}
// generate output file name
string base = b.gen_base_path();
string path = root + "/" + base;
SGPath sgp( path );
sgp.append( "dummy" );
sgp.create_dir( 0755 );
string array_file = path + "/" + b.gen_index_str() + ".arr.gz";
cout << "array_file = " << array_file << endl;
// write the file
gzFile fp;
if ( (fp = gzopen( array_file.c_str(), "wb9" )) == NULL ) {
cout << "ERROR: cannot open " << array_file << " for writing!" << endl;
exit(-1);
}
gzprintf( fp, "%d %d\n", (int)min_x, (int)min_y );
gzprintf( fp, "%d %d %d %d\n", span_x + 1, (int)col_step,
span_y + 1, (int)row_step );
for ( int i = start_x; i <= start_x + span_x; ++i ) {
for ( int j = start_y; j <= start_y + span_y; ++j ) {
gzprintf( fp, "%d ", (int)data[i][j] );
}
gzprintf( fp, "\n" );
}
gzclose(fp);
return true;
}
// write the entire area out in a simple ascii format
bool TGHgt::write_whole_ascii( const string& file ) {
cout << "writing to " << file << endl;
// write the file
gzFile fp;
if ( (fp = gzopen( file.c_str(), "wb9" )) == NULL ) {
cout << "ERROR: cannot open " << file << " for writing!" << endl;
exit(-1);
}
gzprintf( fp, "%d\n%d\n", rows, cols );
for ( int row = rows - 1; row >= 0; row-- ) {
for ( int col = 0; col < cols; col++ ) {
gzprintf( fp, "%d\n", (int)data[col][row] );
}
}
gzclose(fp);
return true;
}
TGHgt::~TGHgt() {
// printf("class TGHgt DEstructor called.\n");
// printf("class TGSrtmBase DEstructor called.\n");
delete [] data;
delete [] output_data;
}
bool
TGHgt::has_non_zero_elev (int start_x, int span_x,
int start_y, int span_y) const
{
for ( int row = start_y; row < start_y + span_y; row++ ) {
for ( int col = start_x; col < start_x + span_x; col++ ) {
if ( data[col][row] != 0 )
return true;
}
}
return false;
}

View file

@ -30,6 +30,8 @@
#include <simgear/compiler.h>
#include "srtmbase.hxx"
#include <zlib.h>
#include <string>
@ -40,31 +42,19 @@
#define MAX_HGT_SIZE 3601
class TGHgt {
class TGHgt : public TGSrtmBase {
private:
// file pointer for input
gzFile fd;
// coordinates (in arc seconds) of south west corner
double originx, originy;
// number of columns and rows
int cols, rows;
// Distance between column and row data points (in arc seconds)
double col_step, row_step;
int hgt_resolution;
// pointers to the actual grid data allocated here
short int (*data)[MAX_HGT_SIZE];
short int (*output_data)[MAX_HGT_SIZE];
int hgt_resolution;
bool remove_tmp_file;
SGPath remove_file_name;
public:
// Constructor, _res must be either "1" for the 1arcsec data or
@ -84,27 +74,7 @@ public:
// load an hgt file
bool load();
// write out the area of data covered by the specified bucket.
// Data is written out column by column starting at the lower left
// hand corner.
bool write_area( const std::string& root, SGBucket& b );
// write the entire area out in a simple ascii format
bool write_whole_ascii( const std::string& file );
// Informational methods
inline double get_originx() const { return originx; }
inline double get_originy() const { return originy; }
inline int get_cols() const { return cols; }
inline int get_rows() const { return rows; }
inline double get_col_step() const { return col_step; }
inline double get_row_step() const { return row_step; }
/**
* Test whether an area contains any non-zero elevations.
*/
bool has_non_zero_elev (int start_x, int span_x,
int start_y, int span_y) const;
virtual short height( int x, int y ) const { return data[x][y]; }
};

151
src/Lib/HGT/srtmbase.cxx Normal file
View file

@ -0,0 +1,151 @@
// hgt.hxx -- SRTM "hgt" data management class
//
// Written by Curtis Olson, started February 2003.
//
// Copyright (C) 2003 Curtis L. Olson - http://www.flightgear.org/~curt
//
// This program is free software; you can redistribute it and/or
// modify it under the terms of the GNU General Public License as
// published by the Free Software Foundation; either version 2 of the
// License, or (at your option) any later version.
//
// This program is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program; if not, write to the Free Software
// Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
//
// $Id: hgt.cxx,v 1.7 2005-12-19 16:06:45 curt Exp $
#ifdef HAVE_CONFIG_H
# include <config.h>
#endif
#include <simgear/compiler.h>
#include <iostream>
#include <zlib.h>
#include "srtmbase.hxx"
using std::cout;
using std::endl;
// write out the area of data covered by the specified bucket. Data
// is written out column by column starting at the lower left hand
// corner.
bool
TGSrtmBase::write_area( const string& root, SGBucket& b ) {
// calculate some boundaries
double min_x = ( b.get_center_lon() - 0.5 * b.get_width() ) * 3600.0;
double max_x = ( b.get_center_lon() + 0.5 * b.get_width() ) * 3600.0;
double min_y = ( b.get_center_lat() - 0.5 * b.get_height() ) * 3600.0;
double max_y = ( b.get_center_lat() + 0.5 * b.get_height() ) * 3600.0;
cout << b << endl;
cout << "width = " << b.get_width() << " height = " << b.get_height()
<< endl;
cout << "min = " << min_x << "," << min_y
<< " max = " << max_x << "," << max_y << endl;
int start_x = (int)((min_x - originx) / col_step);
int span_x = (int)(b.get_width() * 3600.0 / col_step);
int start_y = (int)((min_y - originy) / row_step);
int span_y = (int)(b.get_height() * 3600.0 / row_step);
cout << "start_x = " << start_x << " span_x = " << span_x << endl;
cout << "start_y = " << start_y << " span_y = " << span_y << endl;
// Do a simple sanity checking. But, please, please be nice to
// this write_area() routine and feed it buckets that coincide
// well with the underlying grid structure and spacing.
if ( ( min_x < originx )
|| ( max_x > originx + cols * col_step )
|| ( min_y < originy )
|| ( max_y > originy + rows * row_step ) ) {
cout << " ERROR: bucket at least partially outside HGT data range!" <<
endl;
return false;
}
// If the area is all ocean, skip it.
if ( !has_non_zero_elev(start_x, span_x, start_y, span_y) ) {
cout << "Tile is all zero elevation: skipping" << endl;
return false;
}
// generate output file name
string base = b.gen_base_path();
string path = root + "/" + base;
SGPath sgp( path );
sgp.append( "dummy" );
sgp.create_dir( 0755 );
string array_file = path + "/" + b.gen_index_str() + ".arr.gz";
cout << "array_file = " << array_file << endl;
// write the file
gzFile fp;
if ( (fp = gzopen( array_file.c_str(), "wb9" )) == NULL ) {
cout << "ERROR: cannot open " << array_file << " for writing!" << endl;
exit(-1);
}
gzprintf( fp, "%d %d\n", (int)min_x, (int)min_y );
gzprintf( fp, "%d %d %d %d\n", span_x + 1, (int)col_step,
span_y + 1, (int)row_step );
for ( int i = start_x; i <= start_x + span_x; ++i ) {
for ( int j = start_y; j <= start_y + span_y; ++j ) {
gzprintf( fp, "%d ", (int)height(i,j) );
}
gzprintf( fp, "\n" );
}
gzclose(fp);
return true;
}
// write the entire area out in a simple ascii format
bool TGSrtmBase::write_whole_ascii( const string& file ) {
cout << "writing to " << file << endl;
// write the file
gzFile fp;
if ( (fp = gzopen( file.c_str(), "wb9" )) == NULL ) {
cout << "ERROR: cannot open " << file << " for writing!" << endl;
exit(-1);
}
gzprintf( fp, "%d\n%d\n", rows, cols );
for ( int row = rows - 1; row >= 0; row-- ) {
for ( int col = 0; col < cols; col++ ) {
gzprintf( fp, "%d\n", (int)height(col,row) );
}
}
gzclose(fp);
return true;
}
bool
TGSrtmBase::has_non_zero_elev (int start_x, int span_x,
int start_y, int span_y) const
{
for ( int row = start_y; row < start_y + span_y; row++ ) {
for ( int col = start_x; col < start_x + span_x; col++ ) {
if ( height(col,row) != 0 )
return true;
}
}
return false;
}

82
src/Lib/HGT/srtmbase.hxx Normal file
View file

@ -0,0 +1,82 @@
// hgt.hxx -- SRTM "hgt" data management class
//
// Written by Curtis Olson, started February 2003.
//
// Copyright (C) 2003 Curtis L. Olson - http://www.flightgear.org/~curt
//
// This program is free software; you can redistribute it and/or
// modify it under the terms of the GNU General Public License as
// published by the Free Software Foundation; either version 2 of the
// License, or (at your option) any later version.
//
// This program is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program; if not, write to the Free Software
// Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
//
// $Id: hgt.hxx,v 1.4 2004-11-19 22:25:50 curt Exp $
#ifndef _SRTMBASE_HXX
#define _SRTMBASE_HXX
#ifdef HAVE_CONFIG_H
# include <config.h>
#endif
#include <simgear/compiler.h>
#include <simgear/bucket/newbucket.hxx>
#include <simgear/misc/sg_path.hxx>
class TGSrtmBase {
protected:
// coordinates (in arc seconds) of south west corner
double originx, originy;
// number of columns and rows
int cols, rows;
// Distance between column and row data points (in arc seconds)
double col_step, row_step;
bool remove_tmp_file;
SGPath remove_file_name;
public:
// write out the area of data covered by the specified bucket.
// Data is written out column by column starting at the lower left
// hand corner.
bool write_area( const std::string& root, SGBucket& b );
// write the entire area out in a simple ascii format
bool write_whole_ascii( const std::string& file );
// Informational methods
inline double get_originx() const { return originx; }
inline double get_originy() const { return originy; }
inline int get_cols() const { return cols; }
inline int get_rows() const { return rows; }
inline double get_col_step() const { return col_step; }
inline double get_row_step() const { return row_step; }
/**
* Test whether an area contains any non-zero elevations.
*/
bool has_non_zero_elev (int start_x, int span_x,
int start_y, int span_y) const;
virtual short height( int x, int ) const = 0;
};
#endif // _SRTMBASE_HXX

View file

@ -23,7 +23,7 @@
#---------------------------------------------------------------------------
bin_PROGRAMS = demchop hgtchop fillvoids testassem
bin_PROGRAMS = demchop hgtchop srtmchop fillvoids testassem
demchop_SOURCES = \
demchop.cxx point2d.hxx
@ -41,6 +41,14 @@ hgtchop_LDADD = \
-lsgbucket -lsgmisc -lsgdebug -lsgxml -lz
$(base_LIBS)
srtmchop_SOURCES = \
srtmchop.cxx point2d.hxx
srtmchop_LDADD = \
$(top_builddir)/src/Lib/HGT/libHGT.a \
-lsgbucket -lsgmisc -lsgdebug -lsgxml -lz
$(base_LIBS)
fillvoids_SOURCES = \
fillvoids.cxx

View file

@ -0,0 +1,420 @@
// srtmchop.cxx -- chop up a srtm tiff file into it's corresponding pieces and stuff
// them into the workspace directory
//
// Written by Frederic Bouvier, started February 2009.
//
// Copyright (C) 1997 Curtis L. Olson - http://www.flightgear.org/~curt
//
// This program is free software; you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation; either version 2 of the License, or
// (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program; if not, write to the Free Software
// Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
//
// $Id:$
#ifdef HAVE_CONFIG_H
# include <config.h>
#endif
#include <simgear/compiler.h>
#include <string>
#include <iostream>
#include <iomanip>
#include <fstream>
#include <sstream>
#ifdef _MSC_VER
# include <direct.h>
#endif
#include <simgear/bucket/newbucket.hxx>
#include <simgear/debug/logstream.hxx>
#include <simgear/misc/sg_path.hxx>
#include <tiffio.h>
#include <Polygon/point2d.hxx>
#include <plib/ul.h>
#include <zlib.h>
#include <Lib/HGT/srtmbase.hxx>
using std::cout;
using std::endl;
using std::setfill;
using std::setw;
using std::string;
using std::ifstream;
using std::ostringstream;
using std::ios;
#define MAX_HGT_SIZE 6001
class TGSrtmTiff : public TGSrtmBase {
public:
TGSrtmTiff( const SGPath &file );
~TGSrtmTiff();
bool open( const SGPath &f );
bool close();
// load an hgt file
bool load();
bool is_opened() const { return opened; }
virtual short height( int x, int y ) const { return data[x][y]; }
private:
enum LoadKind { BottomLeft, BottomRight, TopLeft, TopRight };
TGSrtmTiff( const SGPath &file, LoadKind lk );
bool pos_from_name( string name, string &pfx, int &x, int &y );
TIFF* tif;
LoadKind lkind;
string prefix, ext;
SGPath dir;
bool opened;
// pointers to the actual grid data allocated here
short int (*data)[MAX_HGT_SIZE];
short int (*output_data)[MAX_HGT_SIZE];
};
TGSrtmTiff::TGSrtmTiff( const SGPath &file ) {
lkind = BottomLeft;
tif = 0;
remove_tmp_file = false;
data = new short int[MAX_HGT_SIZE][MAX_HGT_SIZE];
output_data = new short int[MAX_HGT_SIZE][MAX_HGT_SIZE];
opened = TGSrtmTiff::open( file );
}
TGSrtmTiff::TGSrtmTiff( const SGPath &file, LoadKind lk ) {
lkind = lk;
tif = 0;
remove_tmp_file = false;
output_data = 0;
if ( lkind == BottomLeft ) {
data = new short int[MAX_HGT_SIZE][MAX_HGT_SIZE];
output_data = new short int[MAX_HGT_SIZE][MAX_HGT_SIZE];
} else if ( lkind == TopLeft ) {
data = new short int[MAX_HGT_SIZE][MAX_HGT_SIZE];
} else if ( lkind == BottomRight ) {
data = new short int[1][MAX_HGT_SIZE];
} else /* if ( lkind == TopRight ) */ {
data = new short int[1][MAX_HGT_SIZE];
}
TGSrtmTiff::open( file );
}
TGSrtmTiff::~TGSrtmTiff() {
delete[] data;
delete[] output_data;
if ( remove_tmp_file ) {
ulDir *dir = ulOpenDir( remove_file_name.dir().c_str() );
if ( dir ) {
ulDirEnt *de;
while ( ( de = ulReadDir( dir ) ) != 0 ) {
if ( !strcmp(de->d_name,".") || !strcmp(de->d_name,"..") || de->d_isdir ) {
continue;
}
SGPath file( remove_file_name.dir() );
file.append( de->d_name );
unlink( file.c_str() );
}
ulCloseDir( dir );
}
rmdir( remove_file_name.dir().c_str() );
}
if ( tif )
TIFFClose( tif );
}
bool TGSrtmTiff::pos_from_name( string name, string &pfx, int &x, int &y ) {
size_t p = name.find( '_' );
if ( p == string::npos )
return false;
pfx = name.substr( 0, p );
name.erase( 0, p + 1 );
p = name.find( '.' );
if ( p == string::npos )
return false;
name.erase( p );
p = name.find( '_' );
if ( p == string::npos )
return false;
x = atoi( name.substr( 0, p ).c_str() );
y = atoi( name.substr( p+1 ).c_str() );
return true;
}
bool TGSrtmTiff::open( const SGPath &f ) {
SGPath file_name = f;
ext = file_name.extension();
dir = file_name.dir();
int x, y;
pos_from_name( file_name.file(), prefix, x, y );
if ( ext == "zip" ) {
// extract the .zip file to /tmp and point the file name
// to the extracted file
SGPath tmp_dir = tempnam( 0, "hgt" );
tmp_dir.append( "dummy" );
tmp_dir.create_dir( 0700 );
cout << "Extracting " << file_name.str() << " to " << tmp_dir.dir() << endl;
string command = "unzip -d \"" + tmp_dir.dir() + "\" " + file_name.str();
system( command.c_str() );
file_name = tmp_dir.dir();
ulDir *dir = ulOpenDir( tmp_dir.dir().c_str() );
if ( dir ) {
ulDirEnt *de;
while ( ( de = ulReadDir( dir ) ) != 0 ) {
if ( !strcmp(de->d_name,".") || !strcmp(de->d_name,"..") || de->d_isdir ) {
continue;
}
SGPath file( de->d_name );
string ext = file.extension();
if ( ext == "TIF" || ext == "tif" ) {
file_name.append( de->d_name );
break;
}
}
ulCloseDir( dir );
}
remove_tmp_file = true;
remove_file_name = file_name.str();
cout << "Proceeding with " << file_name.str() << endl;
}
tif = TIFFOpen( file_name.c_str(), "r" );
if ( !tif ) {
cout << "ERROR: opening " << file_name.str() << " for reading!" << endl;
return false;
}
// Determine originx/originy from file name
originx = ( double( x ) - 37.0 ) * 18000.0;
originy = ( 12.0 - double( y ) ) * 18000.0;
cout << " Origin = " << originx << ", " << originy << endl;
return true;
}
bool TGSrtmTiff::load() {
int size;
cols = rows = size = 6000;
col_step = row_step = 3;
uint32 w, h, d;
uint16 dataType;
uint16 samplesperpixel;
uint16 bitspersample;
TIFFGetField( tif, TIFFTAG_IMAGEWIDTH, &w );
TIFFGetField( tif, TIFFTAG_IMAGELENGTH, &h );
TIFFGetField( tif, TIFFTAG_IMAGEDEPTH, &d );
TIFFGetField( tif, TIFFTAG_SAMPLESPERPIXEL, &samplesperpixel );
TIFFGetField( tif, TIFFTAG_BITSPERSAMPLE, &bitspersample );
TIFFGetField( tif, TIFFTAG_DATATYPE, &dataType );
tdata_t buf = _TIFFmalloc( TIFFScanlineSize( tif ) );
if ( lkind == BottomLeft ) {
uint32 row = 0;
for ( ; row < h; row++ ) {
TIFFReadScanline( tif, buf, row );
uint32 col = 0;
for ( ; col < w; col++ ) {
int16 v = ((int16*)buf)[col];
if ( v == -32768 )
v = 0;
data[col][6000-1-row] = v;
}
for ( ; col < 6000; col++ ) {
data[col][6000-1-row] = 0;
}
}
for ( ; row < 6000; row++ ) {
uint32 col = 0;
for ( ; col < 6000; col++ ) {
data[col][6000-1-row] = 0;
}
}
int x1 = int( originx / 18000.0 ) + 37,
y1 = int( 12 - ( originy / 18000.0 ) ),
x2 = x1 + 1,
y2 = y1 - 1;
if ( x2 > 72 )
x2 -= 72;
{
ostringstream name;
name << prefix << "_" << std::setfill( '0' ) << std::setw( 2 ) << x2 << "_" << std::setfill( '0' ) << std::setw( 2 ) << y1 << "." << ext;
SGPath f = dir;
f.append( name.str() );
if ( f.exists() ) {
TGSrtmTiff s( f.str(), BottomRight );
s.load();
s.close();
for ( int i = 0; i < 6000; ++i ) {
data[6000][i] = s.data[0][i];
}
} else {
for ( int i = 0; i < 6000; ++i ) {
data[6000][i] = 0;
}
}
}
if ( y2 != 0 ) {
ostringstream name;
name << prefix << "_" << std::setfill( '0' ) << std::setw( 2 ) << x1 << "_" << std::setfill( '0' ) << std::setw( 2 ) << y2 << "." << ext;
SGPath f = dir;
f.append( name.str() );
if ( f.exists() ) {
TGSrtmTiff s( f.str(), TopLeft );
s.load();
s.close();
for ( int i = 0; i < 6000; ++i ) {
data[i][6000] = s.data[i][0];
}
} else {
for ( int i = 0; i < 6000; ++i ) {
data[i][6000] = 0;
}
}
} else {
for ( int i = 0; i < 6000; ++i ) {
data[i][6000] = data[i][6000-1];
}
}
if ( y2 != 0 ) {
ostringstream name;
name << prefix << "_" << std::setfill( '0' ) << std::setw( 2 ) << x2 << "_" << std::setfill( '0' ) << std::setw( 2 ) << y2 << "." << ext;
SGPath f = dir;
f.append( name.str() );
if ( f.exists() ) {
TGSrtmTiff s( f.str(), TopRight );
s.load();
s.close();
data[6000][6000] = s.data[0][0];
} else {
data[6000][6000] = 0;
}
} else {
data[6000][6000] = data[6000][6000-1];
}
} else if ( lkind == TopLeft ) {
TIFFReadScanline( tif, buf, 0 );
uint32 col = 0;
for ( ; col < w; col++ ) {
int16 v = ((int16*)buf)[col];
if ( v == -32768 )
v = 0;
data[col][0] = v;
}
for ( ; col < 6000; col++ ) {
data[col][0] = 0;
}
} else if ( lkind == BottomRight ) {
uint32 row = 0;
for ( ; row < h; row++ ) {
TIFFReadScanline( tif, buf, row );
int16 v = ((int16*)buf)[0];
if ( v == -32768 )
v = 0;
data[0][6000-1-row] = v;
}
for ( ; row < 6000; row++ ) {
data[0][6000-1-row] = 0;
}
} else /* if ( lkind == TopRight ) */ {
if ( h == 6000 ) {
TIFFReadScanline( tif, buf, h-1 );
int16 v = ((int16*)buf)[0];
if ( v == -32768 )
v = 0;
data[0][0] = v;
} else {
data[0][0] = 0;
}
}
_TIFFfree(buf);
return true;
}
bool TGSrtmTiff::close() {
if ( tif )
TIFFClose( tif );
tif = 0;
return true;
}
int main(int argc, char **argv) {
sglog().setLogLevels( SG_ALL, SG_WARN );
if ( argc != 3 ) {
cout << "Usage " << argv[0] << " <hgt_file> <work_dir>"
<< endl;
cout << endl;
exit(-1);
}
string hgt_name = argv[1];
string work_dir = argv[2];
SGPath sgp( work_dir );
sgp.append( "dummy" );
sgp.create_dir( 0755 );
TGSrtmTiff hgt( hgt_name );
hgt.load();
hgt.close();
point2d min, max;
min.x = hgt.get_originx() / 3600.0 + SG_HALF_BUCKET_SPAN;
min.y = hgt.get_originy() / 3600.0 + SG_HALF_BUCKET_SPAN;
SGBucket b_min( min.x, min.y );
max.x = (hgt.get_originx() + hgt.get_cols() * hgt.get_col_step()) / 3600.0
- SG_HALF_BUCKET_SPAN;
max.y = (hgt.get_originy() + hgt.get_rows() * hgt.get_row_step()) / 3600.0
- SG_HALF_BUCKET_SPAN;
SGBucket b_max( max.x, max.y );
if ( b_min == b_max ) {
hgt.write_area( work_dir, b_min );
} else {
SGBucket b_cur;
int dx, dy, i, j;
sgBucketDiff(b_min, b_max, &dx, &dy);
cout << "HGT file spans tile boundaries (ok)" << endl;
cout << " dx = " << dx << " dy = " << dy << endl;
if ( (dx > 50) || (dy > 50) ) {
cout << "somethings really wrong!!!!" << endl;
exit(-1);
}
for ( j = 0; j <= dy; j++ ) {
for ( i = 0; i <= dx; i++ ) {
b_cur = sgBucketOffset(min.x, min.y, i, j);
hgt.write_area( work_dir, b_cur );
}
}
}
return 0;
}