geo/tcSrtmData.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 /**
  21  * @file tcSrtmData.cpp
  22  * @brief Block of raw SRTM elevation data.
  23  */
  24
  25 #include "tcSrtmData.h"
  26
  27 #include <QFile>
  28 #include <QDataStream>
  29
  30 #include <cmath>
  31
  32 /*
  33  * Constructors + destructor
  34  */
  35
  36 /// Load SRTM data.
  37 tcSrtmData::tcSrtmData(int lon, int lat, QFile& file)
  38 : m_lon(lon)
  39 , m_lat(lat)
  40 , m_lines(0)
  41 , m_samples(0)
  42 , m_data(0)
  43 {
  44   qint64 elements = file.size() / sizeof(uint16_t);
  45   double side = sqrt(elements);
  46   m_lines = m_samples = (int)side;
  47   Q_ASSERT(m_lines*m_samples == elements);
  48   m_data = new uint16_t[2*m_lines*m_samples];
  49   Q_ASSERT(0 != m_data);
  50
  51   QDataStream in(&file);
  52   in.setByteOrder(QDataStream::BigEndian);
  53   for (int i = 0; i < m_lines; ++i)
  54   {
  55     for (int j = 0; j < m_samples; ++j)
  56     {
  57       Q_ASSERT(!in.atEnd());
  58       uint16_t sample;
  59       in >> sample;
  60       m_data[indexData(0, i, j)] = sample;
  61       if (0x0000 == (0x8000 & sample))
  62       {
  63         m_data[indexData(1, i, j)] = sample;
  64       }
  65       else
  66       {
  67         m_data[indexData(1, i, j)] = (sample & 0x7FFF)+3000;
  68       }
  69     }
  70   }
  71   /// @todo Handle larger files!
  72   Q_ASSERT(in.atEnd());
  73 }
  74
  75 /// Destructor.
  76 tcSrtmData::~tcSrtmData()
  77 {
  78   delete [] m_data;
  79 }
  80
  81 /*
  82  * Accessors
  83  */
  84
  85 /// Get the longitude in degrees.
  86 int tcSrtmData::lon() const
  87 {
  88   return m_lon;
  89 }
  90
  91 /// Get the latitude in degrees.
  92 int tcSrtmData::lat() const
  93 {
  94   return m_lat;
  95 }
  96
  97 /// Get the number of scan lines.
  98 int tcSrtmData::lines() const
  99 {
 100   return m_lines;
 101 }
 102
 103 /// Get the number of samples in each scan line.
 104 int tcSrtmData::samples() const
 105 {
 106   return m_samples;
 107 }
 108
 109 /// Find the geographical angles between samples.
 110 tcGeo tcSrtmData::sampleResolution() const
 111 {
 112   return tcGeo(M_PI/180 / m_samples,
 113                M_PI/180 / m_lines);
 114 }
 115
 116 /// Get the geographical coordinate of a sample.
 117 tcGeo tcSrtmData::sampleCoordinate(int line, int sample) const
 118 {
 119   return tcGeo(M_PI/180*((double)m_lon + (double)sample/m_samples),
 120                M_PI/180*((double)(m_lat+1) - (double)line /m_lines ));
 121 }
 122
 123 /// Does data exist for a sample?
 124 bool tcSrtmData::sampleExists(int line, int sample) const
 125 {
 126   Q_ASSERT(line >= 0 && line < m_lines);
 127   Q_ASSERT(sample >= 0 && sample < m_samples);
 128   return 0x0000 == (0x8000 & m_data[indexData(0, line, sample)]);
 129 }
 130
 131 /// Get the altitude at a sample (interpolating gaps).
 132 uint16_t tcSrtmData::sampleAltitude(int line, int sample, bool corrected, bool* exists) const
 133 {
 134   Q_ASSERT(line >= 0 && line < m_lines);
 135   Q_ASSERT(sample >= 0 && sample < m_samples);
 136   if (0 != exists)
 137   {
 138     *exists = (0x0000 == (0x8000 & m_data[indexData(corrected ? 1 : 0, line, sample)]));
 139   }
 140   return 0x7FFF & m_data[indexData(corrected ? 1 : 0, line, sample)];
 141 }
 142
 143 /// Get the altitude at a coordinate.
 144 uint16_t tcSrtmData::sampleAltitude(const tcGeo& coord, bool corrected, bool* exists) const
 145 {
 146   double dlon = coord.lon()*180.0/M_PI - m_lon;
 147   double dlat = coord.lat()*180.0/M_PI - m_lat;
 148   Q_ASSERT(dlon >= 0.0 && dlon < 1.0);
 149   Q_ASSERT(dlat >= 0.0 && dlat < 1.0);
 150   int line = m_lines - 1 - (int)floor(dlat*m_lines);
 151   int sample = (int)floor(dlon*m_samples);
 152   Q_ASSERT(line >= 0 && line < m_lines);
 153   Q_ASSERT(sample >= 0 && sample < m_samples);
 154   if (0 != exists)
 155   {
 156     *exists = (0x0000 == (0x8000 & m_data[indexData(corrected ? 1 : 0, line, sample)]));
 157   }
 158   return 0x7FFF & m_data[indexData(corrected ? 1 : 0, line, sample)];
 159 }
 160
 161 /*
 162  * Private functions
 163  */
 164
 165 /// Index into the data array.
 166 unsigned int tcSrtmData::indexData(int corrected, int line, int sample) const
 167 {
 168   return (corrected*m_lines + line)*m_samples + sample;
 169 }