arctohgt/arcToHgt.cpp

   1 /***************************************************************************
   2  *   This file is part of Tecorrec.                                        *
   3  *   Copyright 2008 James Hogan <james@albanarts.com>                      *
   4  *                                                                         *
   5  *   Tecorrec is free software: you can redistribute it and/or modify      *
   6  *   it under the terms of the GNU General Public License as published by  *
   7  *   the Free Software Foundation, either version 2 of the License, or     *
   8  *   (at your option) any later version.                                   *
   9  *                                                                         *
  10  *   Tecorrec is distributed in the hope that it will be useful,           *
  11  *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
  12  *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
  13  *   GNU General Public License for more details.                          *
  14  *                                                                         *
  15  *   You should have received a copy of the GNU General Public License     *
  16  *   along with Tecorrec.  If not, write to the Free Software Foundation,  *
  17  *   Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.   *
  18  ***************************************************************************/
  19
  20 #include <QString>
  21 #include <QByteArray>
  22 #include <QFile>
  23 #include <QRegExp>
  24 #include <QtDebug>
  25
  26 #include <cmath>
  27 #include <inttypes.h>
  28
  29 double readDouble(QFile& file)
  30 {
  31   QByteArray line = file.readLine();
  32   static QRegExp intLine("^\\s*\\w+\\s+([-+]?\\d+(\\.\\d+)?)\\s*$");
  33   if (intLine.exactMatch(line))
  34   {
  35     return intLine.cap(1).toDouble();
  36   }
  37   else
  38   {
  39     return -1;
  40   }
  41 }
  42
  43 int readInt(QFile& file)
  44 {
  45   QByteArray line = file.readLine();
  46   static QRegExp intLine("^\\s*\\w+\\s+([-+]?\\d+)\\s*$");
  47   if (intLine.exactMatch(line))
  48   {
  49     return intLine.cap(1).toInt();
  50   }
  51   else
  52   {
  53     return -1;
  54   }
  55 }
  56
  57 int arcToHgt(const QString filename, char* voids)
  58 {
  59   qDebug() << "Reading" << filename;
  60
  61   // first open arcinfo file
  62   QFile af(filename);
  63   if (!af.open(QIODevice::ReadOnly))
  64   {
  65     return 1;
  66   }
  67
  68   // and read the header info
  69   int ncols = readInt(af);
  70   int nrows = readInt(af);
  71   double xllcorner = readDouble(af);
  72   double yllcorner = readDouble(af);
  73   double cellsize = readDouble(af);
  74   int nodata_value = readInt(af);
  75
  76   // find corresponding HGT tiles
  77   int lonMin = (int)floor(xllcorner);
  78   int lonMax = (int)floor(xllcorner + cellsize*(ncols-1));
  79   int latMin = (int)floor(yllcorner);
  80   int latMax = (int)floor(yllcorner + cellsize*(nrows-1));
  81   int tilesX = 1+lonMax-lonMin;
  82   int tilesY = 1+latMax-latMin;
  83   QFile* srtm = new QFile[tilesX*tilesY];
  84   QFile* voidfs = 0;
  85   if (0 != voids)
  86   {
  87     voidfs = new QFile[tilesX*tilesY];
  88   }
  89   for (int lon = lonMin; lon <= lonMax; ++lon)
  90   {
  91     for (int lat = latMin; lat <= latMax; ++lat)
  92     {
  93       QString filename = QString("%3%4%5%6.hgt")
  94                                     .arg(lat >= 0 ? 'N' : 'S')
  95                                     .arg(abs(lat), 2, 10, QLatin1Char('0'))
  96                                     .arg(lon >= 0 ? 'E' : 'W')
  97                                     .arg(abs(lon), 3, 10, QLatin1Char('0'));
  98       int index = tilesX*(lat-latMin) + lon-lonMin;
  99       srtm[index].setFileName(filename);
 100       srtm[index].open(QIODevice::ReadWrite);
 101       if (0 != voids)
 102       {
 103         voidfs[index].setFileName(QString("%1/%2").arg(voids).arg(filename));
 104         voidfs[index].open(QIODevice::ReadOnly);
 105       }
 106       // if it doesn't exist, initialize with void
 107       if (srtm[index].atEnd())
 108       {
 109         qDebug() << "  Creating" << filename;
 110         QDataStream stream(&srtm[index]);
 111         stream.setByteOrder(QDataStream::BigEndian);
 112         uint16_t sample = 0x8000;
 113         for (int i = 0; i < 1201*1201; ++i)
 114         {
 115           stream << sample;
 116         }
 117       }
 118       else
 119       {
 120         qDebug() << "  Writing to" << filename;
 121       }
 122     }
 123   }
 124
 125   // overwrite directly with new values
 126
 127   QTextStream in(&af);
 128   // go through the rows
 129   for (int row = nrows-1; row >= 0; --row)
 130   {
 131     int y = (int)floor(0.5 + yllcorner/cellsize - latMin/cellsize + row);
 132     int latDiv = y/1200;
 133     int latMod = y%1200;
 134     for (int col = 0; col < ncols; ++col)
 135     {
 136       int x = (int)floor(0.5 + xllcorner/cellsize - lonMin/cellsize + col);
 137       int lonDiv = x/1200;
 138       int lonMod = x%1200;
 139
 140       int value;
 141       in >> value;
 142       uint16_t sample = value;
 143
 144       bool inX = (lonDiv >= 0 && lonDiv < tilesX);
 145       bool inY = (latDiv >= 0 && latDiv < tilesY);
 146       bool inXm1 = (lonDiv > 0 && lonDiv <= tilesX);
 147       bool inYm1 = (latDiv > 0 && latDiv <= tilesY);
 148 #define FILEPOS(X,Y) (((1200-(Y))*1201 + (X))*2)
 149       if (inX && inY)
 150       {
 151         int index = latDiv*tilesX + lonDiv;
 152         // if void flag is set, adjust sample
 153         if (0 != voidfs && voidfs[index].isOpen())
 154         {
 155           QDataStream stream(&voidfs[index]);
 156           stream.setByteOrder(QDataStream::BigEndian);
 157           voidfs[index].seek(FILEPOS(lonMod, latMod));
 158           uint16_t mask;
 159           stream >> mask;
 160           if (mask & 0x8000)
 161           {
 162             sample |= 0x8000;
 163           }
 164         }
 165         QDataStream stream(&srtm[index]);
 166         stream.setByteOrder(QDataStream::BigEndian);
 167         srtm[index].seek(FILEPOS(lonMod, latMod));
 168         stream << sample;
 169       }
 170       bool edgeX = (lonMod == 0);
 171       bool edgeY = (latMod == 0);
 172       if (edgeX && inXm1 && inY)
 173       {
 174         int index = latDiv*tilesX + lonDiv-1;
 175         QDataStream stream(&srtm[index]);
 176         stream.setByteOrder(QDataStream::BigEndian);
 177         srtm[index].seek(FILEPOS(1200, latMod));
 178         stream << sample;
 179       }
 180       if (edgeY && inX && inYm1)
 181       {
 182         int index = (latDiv-1)*tilesX + lonDiv;
 183         QDataStream stream(&srtm[index]);
 184         stream.setByteOrder(QDataStream::BigEndian);
 185         srtm[index].seek(FILEPOS(lonMod, 1200));
 186         stream << sample;
 187       }
 188       if (edgeX && edgeY && inXm1 && inYm1)
 189       {
 190         int index = (latDiv-1)*tilesX + lonDiv-1;
 191         QDataStream stream(&srtm[index]);
 192         stream.setByteOrder(QDataStream::BigEndian);
 193         srtm[index].seek(FILEPOS(1200, 1200));
 194         stream << sample;
 195       }
 196     }
 197   }
 198
 199   delete [] srtm;
 200   if (0 != voidfs)
 201   {
 202     delete [] voidfs;
 203   }
 204
 205   return 0;
 206 }
 207
 208 int main(int argc, char** argv)
 209 {
 210   qDebug() << "Welcome to arcToHgt";
 211   int error = 0;
 212   // Convert the file in each argument to HGT, overwriting any existing data
 213   char * voids = 0;
 214   for (int i = 1; i < argc; ++i)
 215   {
 216     if (0 == strcmp(argv[i], "--voids"))
 217     {
 218       // next argument is voids directory
 219       if (++i < argc)
 220       {
 221         voids = argv[i];
 222         qDebug() << "Using voids from " << voids;
 223       }
 224     }
 225     else
 226     {
 227       error |= arcToHgt(argv[i], voids);
 228     }
 229   }
 230   return error;
 231 }