source: ntrip/trunk/BNC/src/bncantex.cpp @ 7865

Last change on this file since 7865 was 7865, checked in by stuerze, 18 months ago

minor changes

File size: 11.2 KB
Line 
1// Part of BNC, a utility for retrieving decoding and
2// converting GNSS data streams from NTRIP broadcasters.
3//
4// Copyright (C) 2007
5// German Federal Agency for Cartography and Geodesy (BKG)
6// http://www.bkg.bund.de
7// Czech Technical University Prague, Department of Geodesy
8// http://www.fsv.cvut.cz
9//
10// Email: euref-ip@bkg.bund.de
11//
12// This program is free software; you can redistribute it and/or
13// modify it under the terms of the GNU General Public License
14// as published by the Free Software Foundation, version 2.
15//
16// This program is distributed in the hope that it will be useful,
17// but WITHOUT ANY WARRANTY; without even the implied warranty of
18// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19// GNU General Public License for more details.
20//
21// You should have received a copy of the GNU General Public License
22// along with this program; if not, write to the Free Software
23// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
24
25/* -------------------------------------------------------------------------
26 * BKG NTRIP Client
27 * -------------------------------------------------------------------------
28 *
29 * Class:      bncAntex
30 *
31 * Purpose:    Antenna Phase Centers and Variations from ANTEX File
32 *
33 * Author:     L. Mervart
34 *
35 * Created:    26-Jan-2011
36 *
37 * Changes:
38 *
39 * -----------------------------------------------------------------------*/
40
41#include <iostream>
42#include <newmatio.h>
43
44#include "bncantex.h"
45#include "pppModel.h"
46
47using namespace std;
48
49// Constructor
50////////////////////////////////////////////////////////////////////////////
51bncAntex::bncAntex() {
52}
53
54// Constructor
55////////////////////////////////////////////////////////////////////////////
56bncAntex::bncAntex(const char* fileName) {
57  readFile(QString(fileName));
58}
59
60// Destructor
61////////////////////////////////////////////////////////////////////////////
62bncAntex::~bncAntex() {
63  QMapIterator<QString, t_antMap*> it(_maps);
64  while (it.hasNext()) {
65    it.next();
66    delete it.value();
67  }
68}
69
70// Print
71////////////////////////////////////////////////////////////////////////////
72void bncAntex::print() const {
73  QMapIterator<QString, t_antMap*> itAnt(_maps);
74  while (itAnt.hasNext()) {
75    itAnt.next();
76    t_antMap* map = itAnt.value();
77    cout << map->antName.toAscii().data() << endl;
78    cout << "    " << map->zen1 << " " << map->zen2 << " " << map->dZen << endl;
79    QMapIterator<t_frequency::type, t_frqMap*> itFrq(map->frqMap);
80    while (itFrq.hasNext()) {
81      itFrq.next();
82      const t_frqMap* frqMap = itFrq.value();
83      cout << "    " << frqMap->neu[0] << " "
84                     << frqMap->neu[1] << " "
85                     << frqMap->neu[2] << endl;
86      cout << "    " << frqMap->pattern.t();
87    }
88    cout << endl;
89  }
90}
91
92// Print
93////////////////////////////////////////////////////////////////////////////
94QString bncAntex::pcoSinexString(const std::string& antName, t_frequency::type frqType) {
95
96  if (antName.find("NULLANTENNA") != string::npos) {
97    return QString(" ------ ------ ------");
98  }
99
100  QString antNameQ = antName.c_str();
101  if (_maps.find(antNameQ) == _maps.end()) {
102    return QString(" ------ ------ ------");
103  }
104
105  t_antMap* map = _maps[antNameQ];
106  if (map->frqMap.find(frqType) == map->frqMap.end()) {
107    return QString(" ------ ------ ------");
108  }
109
110  t_frqMap* frqMap = map->frqMap[frqType];
111
112  QString u, n,e;
113  u.sprintf("%+6.4f" ,frqMap->neu[2]); if (u.mid(1,1) == "0") {u.remove(1,1);}
114  n.sprintf("%+6.4f" ,frqMap->neu[0]); if (n.mid(1,1) == "0") {n.remove(1,1);}
115  e.sprintf("%+6.4f" ,frqMap->neu[1]); if (e.mid(1,1) == "0") {e.remove(1,1);}
116
117  return QString(" %1 %2 %3").arg(u).arg(n).arg(e);
118}
119
120// Read ANTEX File
121////////////////////////////////////////////////////////////////////////////
122t_irc bncAntex::readFile(const QString& fileName) {
123
124  QFile inFile(fileName);
125  inFile.open(QIODevice::ReadOnly | QIODevice::Text);
126
127  QTextStream in(&inFile);
128
129  t_antMap* newAntMap = 0;
130  t_frqMap* newFrqMap = 0;
131
132  while ( !in.atEnd() ) {
133    QString line = in.readLine();
134
135    // Start of Antenna
136    // ----------------
137    if      (line.indexOf("START OF ANTENNA") == 60) {
138      if (newAntMap) {
139        delete newAntMap;
140        return failure;
141      }
142      else {
143        delete newAntMap;
144        newAntMap = new t_antMap();
145      }
146    }
147
148    // End of Antenna
149    // --------------
150    else if (line.indexOf("END OF ANTENNA") == 60) {
151      if (newAntMap) {
152        if (_maps.contains(newAntMap->antName)) {
153          delete _maps[newAntMap->antName];
154        }
155        _maps[newAntMap->antName] = newAntMap;
156        newAntMap = 0;
157      }
158      else {
159        delete newAntMap;
160        return failure;
161      }
162    }
163
164    // Antenna Reading in Progress
165    // ---------------------------
166    else if (newAntMap) {
167      if      (line.indexOf("TYPE / SERIAL NO") == 60) {
168        if (line.indexOf("BLOCK I") == 0 ||
169            line.indexOf("GLONASS") == 0 ||
170            line.indexOf("QZSS") == 0 ||
171            line.indexOf("BEIDOU") == 0 ||
172            line.indexOf("GALILEO") == 0 ||
173            line.indexOf("IRNSS") == 0 ){
174          newAntMap->antName = line.mid(20,3);
175        }
176        else {
177          newAntMap->antName = line.mid(0,20);
178        }
179      }
180      else if (line.indexOf("ZEN1 / ZEN2 / DZEN") == 60) {
181        QTextStream inLine(&line, QIODevice::ReadOnly);
182        inLine >> newAntMap->zen1 >> newAntMap->zen2 >> newAntMap->dZen;
183      }
184
185      // Start of Frequency
186      // ------------------
187      else if (line.indexOf("START OF FREQUENCY") == 60) {
188        if (newFrqMap) {
189          delete newFrqMap;
190          delete newAntMap;
191          return failure;
192        }
193        else {
194          newFrqMap = new t_frqMap();
195        }
196      }
197
198      // End of Frequency
199      // ----------------
200      else if (line.indexOf("END OF FREQUENCY") == 60) {
201        if (newFrqMap) {
202          t_frequency::type frqType = t_frequency::dummy;
203          if      (line.indexOf("G01") == 3) {
204            frqType = t_frequency::G1;
205          }
206          else if (line.indexOf("G02") == 3) {
207            frqType = t_frequency::G2;
208          }
209          else if (line.indexOf("R01") == 3) {
210            frqType = t_frequency::R1;
211          }
212          else if (line.indexOf("R02") == 3) {
213            frqType = t_frequency::R2;
214          }
215          else if (line.indexOf("E01") == 3) {
216            frqType = t_frequency::E1;
217          }
218          else if (line.indexOf("E05") == 3) {
219            frqType = t_frequency::E5;
220          }
221          else if (line.indexOf("E06") == 3) {
222            frqType = t_frequency::E6;
223          }
224          else if (line.indexOf("E07") == 3) {
225            frqType = t_frequency::E7;
226          }
227          else if (line.indexOf("E08") == 3) {
228            frqType = t_frequency::E8;
229          }
230          else if (line.indexOf("J01") == 3) {
231            frqType = t_frequency::J1;
232          }
233          else if (line.indexOf("J02") == 3) {
234            frqType = t_frequency::J2;
235          }
236          else if (line.indexOf("J05") == 3) {
237            frqType = t_frequency::J5;
238          }
239          else if (line.indexOf("J06") == 3) {
240            frqType = t_frequency::J6;
241          }
242          else if (line.indexOf("C02") == 3) {
243            frqType = t_frequency::C2;
244          }
245          else if (line.indexOf("C06") == 3) {
246            frqType = t_frequency::C6;
247          }
248          else if (line.indexOf("C07") == 3) {
249            frqType = t_frequency::C7;
250          }
251          if (frqType != t_frequency::dummy) {
252            if (newAntMap->frqMap.find(frqType) != newAntMap->frqMap.end()) {
253              delete newAntMap->frqMap[frqType];
254            }
255            newAntMap->frqMap[frqType] = newFrqMap;
256          }
257          else {
258            delete newFrqMap;
259          }
260          newFrqMap = 0;
261        }
262        else {
263          delete newAntMap;
264          return failure;
265        }
266      }
267
268      // Frequency Reading in Progress
269      // -----------------------------
270      else if (newFrqMap) {
271        if      (line.indexOf("NORTH / EAST / UP") == 60) {
272          QTextStream inLine(&line, QIODevice::ReadOnly);
273          inLine >> newFrqMap->neu[0] >> newFrqMap->neu[1] >> newFrqMap->neu[2];
274          newFrqMap->neu[0] *= 1e-3;
275          newFrqMap->neu[1] *= 1e-3;
276          newFrqMap->neu[2] *= 1e-3;
277        }
278        else if (line.indexOf("NOAZI") == 3) {
279          QTextStream inLine(&line, QIODevice::ReadOnly);
280          int nPat = int((newAntMap->zen2-newAntMap->zen1)/newAntMap->dZen) + 1;
281          newFrqMap->pattern.ReSize(nPat);
282          QString dummy;
283          inLine >> dummy;
284          for (int ii = 0; ii < nPat; ii++) {
285            inLine >> newFrqMap->pattern[ii];
286          }
287          newFrqMap->pattern *= 1e-3;
288        }
289      }
290    }
291  }
292  inFile.close();
293  delete newFrqMap;
294  delete newAntMap;
295
296  return success;
297}
298
299// Satellite Antenna Offset
300////////////////////////////////////////////////////////////////////////////
301t_irc bncAntex::satCoMcorrection(const QString& prn, double Mjd,
302                                 const ColumnVector& xSat, ColumnVector& dx) {
303
304  t_frequency::type frqType = t_frequency::dummy;
305
306  if      (prn[0] == 'G') {
307    frqType = t_frequency::G1;
308  }
309  else if (prn[0] == 'R') {
310    frqType = t_frequency::R1;
311  }
312
313  QMap<QString, t_antMap*>::const_iterator it = _maps.find(prn.mid(0,3));
314  if (it != _maps.end()) {
315    t_antMap* map = it.value();
316    if (map->frqMap.find(frqType) != map->frqMap.end()) {
317
318      double* neu = map->frqMap[frqType]->neu;
319
320      // Unit Vectors sz, sy, sx
321      // -----------------------
322      ColumnVector sz = -xSat;
323      sz /= sqrt(DotProduct(sz,sz));
324
325      ColumnVector xSun = BNC_PPP::t_astro::Sun(Mjd);
326      xSun /= sqrt(DotProduct(xSun,xSun));
327
328      ColumnVector sy = crossproduct(sz, xSun);
329      sy /= sqrt(DotProduct(sy,sy));
330
331      ColumnVector sx = crossproduct(sy, sz);
332
333      dx[0] = sx[0] * neu[0] + sy[0] * neu[1] + sz[0] * neu[2];
334      dx[1] = sx[1] * neu[0] + sy[1] * neu[1] + sz[1] * neu[2];
335      dx[2] = sx[2] * neu[0] + sy[2] * neu[1] + sz[2] * neu[2];
336
337      return success;
338    }
339  }
340
341  return failure;
342}
343
344//
345////////////////////////////////////////////////////////////////////////////
346double bncAntex::rcvCorr(const string& antName, t_frequency::type frqType,
347                         double eleSat, double azSat, bool& found) const {
348
349  if (antName.find("NULLANTENNA") != string::npos) {
350    found = true;
351    return 0.0;
352  }
353
354  QString antNameQ = antName.c_str();
355
356  if (_maps.find(antNameQ) == _maps.end()) {
357    found = false;
358    return 0.0;
359  }
360
361  t_antMap* map = _maps[antNameQ];
362  if (map->frqMap.find(frqType) == map->frqMap.end()) {
363    found = false;
364    return 0.0;
365  }
366
367  t_frqMap* frqMap = map->frqMap[frqType];
368
369  double var = 0.0;
370  if (frqMap->pattern.ncols() > 0) {
371    double zenDiff = 999.999;
372    double zenSat  = 90.0 - eleSat * 180.0 / M_PI;
373    unsigned iZen = 0;
374    for (double zen = map->zen1; zen <= map->zen2; zen += map->dZen) {
375      iZen += 1;
376      double newZenDiff = fabs(zen - zenSat);
377      if (newZenDiff < zenDiff) {
378        zenDiff = newZenDiff;
379        var = frqMap->pattern(iZen);
380      }
381    }
382  }
383
384  found = true;
385  return var - frqMap->neu[0] * cos(azSat)*cos(eleSat)
386             - frqMap->neu[1] * sin(azSat)*cos(eleSat)
387             - frqMap->neu[2] * sin(eleSat);
388
389}
Note: See TracBrowser for help on using the repository browser.