From 02d9174161f2c6e0244e873e9cc5378252da4588 Mon Sep 17 00:00:00 2001 From: James Hogan Date: Wed, 4 Mar 2009 22:27:01 +0000 Subject: [PATCH] Simple program to read arcinfo ascii files and convert into hgt elevation data (for cgair-csi dataset) --- CMakeLists.txt | 1 + arctohgt/CMakeLists.txt | 23 +++++++ arctohgt/arcToHgt.cpp | 180 ++++++++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 204 insertions(+) create mode 100644 arctohgt/CMakeLists.txt create mode 100644 arctohgt/arcToHgt.cpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 029a720..e75fbf8 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -16,6 +16,7 @@ include_directories(${QT_INCLUDES} ) add_subdirectory(maths) add_subdirectory(geo) +add_subdirectory(arctohgt) set(tecorrec_SRCS main.cpp diff --git a/arctohgt/CMakeLists.txt b/arctohgt/CMakeLists.txt new file mode 100644 index 0000000..b16dca1 --- /dev/null +++ b/arctohgt/CMakeLists.txt @@ -0,0 +1,23 @@ +cmake_minimum_required(VERSION 2.6) +project(arctohgt) + +find_package(Qt4 REQUIRED) + +include(${QT_USE_FILE}) +include_directories(${QT_INCLUDES}) + +set(arctohgt_SRCS + arcToHgt.cpp +) + +set(arctohgt_HEADERS +) + +#qt4_wrap_cpp(arctohgt_MOCS ${arctohgt_HEADERS}) + +add_executable(arctohgt ${arctohgt_SRCS} ${arctohgt_MOCS}) + +target_link_libraries(arctohgt + ${QT_LIBRARIES} +) + diff --git a/arctohgt/arcToHgt.cpp b/arctohgt/arcToHgt.cpp new file mode 100644 index 0000000..bec96e8 --- /dev/null +++ b/arctohgt/arcToHgt.cpp @@ -0,0 +1,180 @@ +/*************************************************************************** + * This file is part of Tecorrec. * + * Copyright 2008 James Hogan * + * * + * Tecorrec 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. * + * * + * Tecorrec 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 Tecorrec. If not, write to the Free Software Foundation, * + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. * + ***************************************************************************/ + +#include +#include +#include +#include +#include + +#include + +double readDouble(QFile& file) +{ + QByteArray line = file.readLine(); + static QRegExp intLine("^\\s*\\w+\\s+([-+]?\\d+(\\.\\d+)?)\\s*$"); + if (intLine.exactMatch(line)) + { + return intLine.cap(1).toDouble(); + } + else + { + return -1; + } +} + +int readInt(QFile& file) +{ + QByteArray line = file.readLine(); + static QRegExp intLine("^\\s*\\w+\\s+([-+]?\\d+)\\s*$"); + if (intLine.exactMatch(line)) + { + return intLine.cap(1).toInt(); + } + else + { + return -1; + } +} + +int arcToHgt(const QString filename) +{ + qDebug() << "Reading" << filename; + + // first open arcinfo file + QFile af(filename); + af.open(QIODevice::ReadOnly); + + // and read the header info + int ncols = readInt(af); + int nrows = readInt(af); + double xllcorner = readDouble(af); + double yllcorner = readDouble(af); + double cellsize = readDouble(af); + int nodata_value = readInt(af); + + // find corresponding HGT tiles + int lonMin = (int)floor(xllcorner); + int lonMax = (int)floor(xllcorner + cellsize*(ncols-1)); + int latMin = (int)floor(yllcorner); + int latMax = (int)floor(yllcorner + cellsize*(nrows-1)); + int tilesX = 1+lonMax-lonMin; + int tilesY = 1+latMax-latMin; + QFile* srtm = new QFile[tilesX*tilesY]; + for (int lon = lonMin; lon <= lonMax; ++lon) + { + for (int lat = latMin; lat <= latMax; ++lat) + { + QString filename = QString("%3%4%5%6.hgt") + .arg(lat >= 0 ? 'N' : 'S') + .arg(abs(lat), 2, 10, QLatin1Char('0')) + .arg(lon >= 0 ? 'E' : 'W') + .arg(abs(lon), 3, 10, QLatin1Char('0')); + int index = tilesX*(lat-latMin) + lon-lonMin; + srtm[index].setFileName(filename); + srtm[index].open(QIODevice::WriteOnly); + // if it doesn't exist, initialize with void + if (srtm[index].atEnd()) + { + QDataStream stream(&srtm[index]); + stream.setByteOrder(QDataStream::BigEndian); + uint16_t sample = 0x8000; + for (int i = 0; i < 1201*1201; ++i) + { + stream << sample; + } + } + } + } + + // overwrite directly with new values + + QTextStream in(&af); + // go through the rows + for (int row = 0; row < nrows; ++row) + { + int y = (int)floor(0.5 + yllcorner/cellsize - latMin/cellsize + nrows-1-row); + int latDiv = y/1200; + int latMod = y%1200; + for (int col = 0; col < ncols; ++col) + { + int x = (int)floor(0.5 + xllcorner/cellsize - lonMin/cellsize + col); + int lonDiv = x/1200; + int lonMod = x%1200; + + int value; + in >> value; + uint16_t sample = value; + + bool inX = (lonDiv >= 0 && lonDiv < tilesX); + bool inY = (latDiv >= 0 && latDiv < tilesY); + bool inXm1 = (lonDiv > 0 && lonDiv <= tilesX); + bool inYm1 = (latDiv > 0 && latDiv <= tilesY); +#define FILEPOS(X,Y) (((1201-(Y))*1201 + (X))*2) + if (inX && inY) + { + int index = latDiv*tilesX + lonDiv; + QDataStream stream(&srtm[index]); + stream.setByteOrder(QDataStream::BigEndian); + srtm[index].seek(FILEPOS(lonMod, latMod)); + stream << sample; + } + bool edgeX = (lonMod == 0); + bool edgeY = (latMod == 0); + if (edgeX && inXm1 && inY) + { + int index = latDiv*tilesX + lonDiv-1; + QDataStream stream(&srtm[index]); + stream.setByteOrder(QDataStream::BigEndian); + srtm[index].seek(FILEPOS(1200, latMod)); + stream << sample; + } + if (edgeY && inX && inYm1) + { + int index = (latDiv-1)*tilesX + lonDiv; + QDataStream stream(&srtm[index]); + stream.setByteOrder(QDataStream::BigEndian); + srtm[index].seek(FILEPOS(lonMod, 1200)); + stream << sample; + } + if (edgeX && edgeY && inXm1 && inYm1) + { + int index = (latDiv-1)*tilesX + lonDiv-1; + QDataStream stream(&srtm[index]); + stream.setByteOrder(QDataStream::BigEndian); + srtm[index].seek(FILEPOS(1200, 1200)); + stream << sample; + } + } + } + + delete [] srtm; +} + +int main(int argc, char** argv) +{ + qDebug() << "Welcome to arcToHgt"; + int error = 0; + // Convert the file in each argument to HGT, overwriting any existing data + for (int i = 1; i < argc; ++i) + { + error |= arcToHgt(argv[i]); + } + return error; +} -- 2.11.4.GIT