C++ - EOP（地球姿勢パラメータ）データファイル 生成！

IERS（International Earth Rotation and Reference systems Service; 国際地球回転観測事業） の EOP（Earth Orientation Parameter; 地球姿勢（回転）パラメータ）から確定／速報／推定値を抽出し、テキストファイルを生成するスクリプトを C++ で作成しました。（人工衛星の位置を正確に計算するために使用するデータ）

過去には Ruby や Python で生成しました。（但し、 CSV ファイル）

• Ruby, Python - EOP（地球姿勢パラメータ）CSV 生成！

0. 前提条件

• Debian GNU/Linux 10.8 (64bit) での作業を想定。
• GCC 10.2.0 (G++ 10.2.0) (C++17) でのコンパイルを想定。
• ここでは EOP（Earth Orientation Parameter; 地球姿勢（回転）パラメータ）が何かについての説明はしない。

1. 事前準備

今回使用するデータを用意しておく。

• こちら から /standard/finals2000A.all, /daily/finals2000A.daily をダウンロードし、 以下で作成するソースコードと同じディレクトリ配下に配置する。

2. C++ ソースコードの作成

File: make_eop.rb

/***********************************************************
  IERS の Buttelin A テキストデータから 
  EOP(Polar Motion etc.) テキストファイルを生成
  * 予め [こちら](ftp://ftp.iers.org/products/eop/rapid/) か
    らダウンロードしておいたものを使用する。
    (IAU 2000A 章動理論によるデータ 
     finals2000A.all", "finals2000A.daily")
  * 1日のデータに速報値(区分"I")と確定値がある場合は、確定値
    を優先。
  * 2ファイルで重複する日付のデータは "finals2000A.daily" を
    優先。
  * 取得する項目：
      date, flag_pm, pm_x, pm_x_e, pm_y, pm_y_e,
      flag_dut, dut1, dut1_e, lod, lod_e,
      flag_nut, nut_x, nut_x_e, nut_y, nut_y_e

    DATE        AUTHOR       VERSION
    2021.05.27  mk-mode.com  1.00 新規作成

  Copyright(C) 2021 mk-mode.com All Rights Reserved.
------------------------------------------------------------
  引数 : なし
------------------------------------------------------------
  $ g++102 -std=c++17 -Wall -O2 --pedantic-errors -o make_eop make_eop.cpp
***********************************************************/
#include <cstdlib>   // for EXIT_XXXX
#include <fstream>
#include <iomanip>
#include <iostream>
#include <map>
#include <string>

namespace make_eop {

// 定数
static constexpr char kFAll[]   = "finals2000A.all";
static constexpr char kFDaily[] = "finals2000A.daily";
static constexpr char kFEop[]   = "eop.txt";
static constexpr char kS6[]     = "      ";
static constexpr char kS7[]     = "       ";
static constexpr char kS8[]     = "        ";
static constexpr char kS9[]     = "         ";
static constexpr char kS10[]    = "          ";
static constexpr char kZ3[]     = "   0.000";
static constexpr char kZ6[]     = "0.000000";
static constexpr char kZ7[]     = "0.0000000";

/*
 * @brief       データ取得
 *
 * @param[in]   ファイル名 (char)
 * @param[ref]  データ (map<string, string>)
 * @return      true|false
 */
bool get_data(std::string f, std::map<std::string, std::string>& data) {
  std::string  l;
  std::string  s;
  char         flag_pm;
  char         flag_dut;
  char         flag_nut;
  unsigned int y;
  unsigned int m;
  unsigned int d;
  std::string  date;
  std::string  mjd;
  std::string  pm_x;
  std::string  pm_x_e;
  std::string  pm_y;
  std::string  pm_y_e;
  std::string  pm_x_f;
  std::string  pm_y_f;
  std::string  dut1;
  std::string  dut1_e;
  std::string  lod;
  std::string  lod_e;
  std::string  dut1_f;
  std::string  nut_x;
  std::string  nut_x_e;
  std::string  nut_y;
  std::string  nut_y_e;
  std::string  nut_x_f;
  std::string  nut_y_f;
  std::stringstream ss;

  try {
    // File open
    std::ifstream ifs(f);
    if (!ifs) return false;

    // File read
    while (getline(ifs, l)) {
      flag_pm  = l[16];
      flag_dut = l[57];
      flag_nut = l[95];
      y        = stoi(l.substr(0, 2));
      m        = stoi(l.substr(2, 2));
      d        = stoi(l.substr(4, 2));
      mjd      = l.substr(  7,  8);
      pm_x     = l.substr( 18,  9);
      pm_x_e   = l.substr( 28,  8);
      pm_y     = l.substr( 37,  9);
      pm_y_e   = l.substr( 47,  8);
      dut1     = l.substr( 58, 10);
      dut1_e   = l.substr( 69,  9);
      lod      = l.substr( 79,  7);
      lod_e    = l.substr( 87,  6);
      nut_x    = l.substr( 97,  9);
      nut_x_e  = l.substr(107,  8);
      nut_y    = l.substr(116,  9);
      nut_y_e  = l.substr(126,  8);
      pm_x_f   = l.substr(135,  9);
      pm_y_f   = l.substr(145,  9);
      dut1_f   = l.substr(155, 10);
      nut_x_f  = l.substr(167,  8);
      nut_y_f  = l.substr(177,  8);

      // Flag check
      if (f == kFAll &&
          flag_pm  == 'P' && flag_dut == 'P' && flag_nut == 'P') {
        break;
      }

      // Date
      if (y > 72) { y += 1900; } else { y += 2000; }
      ss.str("");
      ss.clear();
      ss << std::setfill('0')
         << std::setw(2) << y << "-"
         << std::setw(2) << m << "-"
         << std::setw(2) << d << std::setfill(' ');
      date = ss.str();

      // Polar Motion
      if (pm_x != kS9) {
        ss.str("");
        ss.clear();
        ss << std::setw(9) << std::fixed << std::setprecision(6)
           << stod(pm_x);
        pm_x = ss.str();
      }
      if (pm_x_e != kS8) {
        ss.str("");
        ss.clear();
        ss << std::setw(8) << std::fixed << std::setprecision(6)
           << stod(pm_x_e);
        pm_x_e = ss.str();
      }
      if (pm_y != kS9) {
        ss.str("");
        ss.clear();
        ss << std::setw(9) << std::fixed << std::setprecision(6)
           << stod(pm_y);
        pm_y = ss.str();
      }
      if (pm_y_e != kS8) {
        ss.str("");
        ss.clear();
        ss << std::setw(8) << std::fixed << std::setprecision(6)
           << stod(pm_y_e);
        pm_y_e = ss.str();
      }
      if (pm_x_f != kS9) {
        ss.str("");
        ss.clear();
        ss << std::setw(9) << std::fixed << std::setprecision(6)
           << stod(pm_x_f);
        pm_x_f = ss.str();
      }
      if (pm_y_f != kS9) {
        ss.str("");
        ss.clear();
        ss << std::setw(9) << std::fixed << std::setprecision(6)
           << stod(pm_y_f);
        pm_y_f = ss.str();
      }
      if (pm_x_f != kS9 && pm_y_f != kS9) { flag_pm = 'F'; }
      if (pm_x_f != kS9) {
        pm_x   = pm_x_f;
        pm_x_e = kZ6;
      }
      if (pm_y_f != kS9) {
        pm_y   = pm_y_f;
        pm_y_e = kZ6;
      }

      // UT1 - UTC, LOD
      if (lod != kS7) {
        ss.str("");
        ss.clear();
        ss << std::setw(7) << std::fixed << std::setprecision(4)
           << stod(lod);
        lod = ss.str();
      }
      if (lod_e != kS6) {
        ss.str("");
        ss.clear();
        ss << std::setw(6) << std::fixed << std::setprecision(4)
           << stod(lod_e);
        lod_e = ss.str();
      }
      if (dut1_f != kS10) {
        ss.str("");
        ss << std::setw(10) << std::fixed << std::setprecision(7)
           << stod(dut1_f);
        dut1_f = ss.str();
      }

      if (dut1_f != kS10) {
        flag_dut = 'F';
        dut1     = dut1_f;
        dut1_e   = kZ7;
      }

      // Nutation
      if (nut_x != kS9) {
        ss.str("");
        ss.clear();
        ss << std::setw(9) << std::fixed << std::setprecision(3)
           << stod(nut_x);
        nut_x = ss.str();
      }
      if (nut_x_e != kS8) {
        ss.str("");
        ss.clear();
        ss << std::setw(8) << std::fixed << std::setprecision(3)
           << stod(nut_x_e);
        nut_x_e = ss.str();
      }
      if (nut_y != kS9) {
        ss.str("");
        ss.clear();
        ss << std::setw(9) << std::fixed << std::setprecision(3)
           << stod(nut_y);
        nut_y = ss.str();
      }
      if (nut_y_e != kS8) {
        ss.str("");
        ss.clear();
        ss << std::setw(8) << std::fixed << std::setprecision(3)
           << stod(nut_y_e);
        nut_y_e = ss.str();
      }
      if (nut_x_f != kS8) {
        ss.str("");
        ss.clear();
        ss << std::setw(8) << std::fixed << std::setprecision(3)
           << stod(nut_x_f);
        nut_x_f = ss.str();
      }
      if (nut_y_f != kS8) {
        ss.str("");
        ss.clear();
        ss << std::setw(8) << std::fixed << std::setprecision(3)
           << stod(nut_y_f);
        nut_y_f = ss.str();
      }
      if (nut_x_f != kS8 && nut_y_f != kS8) { flag_nut = 'F'; }
      if (nut_x_f != kS8) {
        nut_x   = " " + nut_x_f;
        nut_x_e = kZ3;
      }
      if (nut_y_f != kS8) {
        nut_y   = " " + nut_y_f;
        nut_y_e = kZ3;
      }

      // 1行整形
      data[date] =
        mjd      + " " +
        flag_pm  + " " +
        pm_x     + " " + pm_x_e  + " " +
        pm_y     + " " + pm_y_e  + " " +
        flag_dut + " " +
        dut1     + " " + dut1_e  + " " +
        lod      + " " + lod_e   + " " +
        flag_nut + " " +
        nut_x    + " " + nut_x_e + " " +
        nut_y    + " " + nut_y_e;
    }
  } catch (...) {
    return false;
  }

  return true;
}

/*
 * @brief       テキストファイル書き込み
 *
 * @param[ref]  データ (map<string, string>)
 * @return      true|false
 */
bool write_data(std::map<std::string, std::string>& data) {
  std::string f(kFEop);

  try {
    // File open
    std::ofstream ofs(f);
    if (!ofs) return 0;

    // 書き込み
    for (auto& [k, v]: data) { ofs << k + " " + v << std::endl; }

    // File close
    ofs.close();
  } catch (...) {
    return false;
  }

  return true;
}

}  // namespace make_eop

int main(int argc, char* argv[]) {
  namespace ns = make_eop;
  std::map<std::string, std::string> data;

  try {
    // "finals2000A.all" 読み込み
    if (!ns::get_data(ns::kFAll, data)) {
      std::cout << "[ERROR] Could not read '" << ns::kFAll << "'!" << std::endl;
      return EXIT_FAILURE;
    }

    // "finals2000A.daily" 読み込み
    if (!ns::get_data(ns::kFDaily, data)) {
      std::cout << "[ERROR] Could not read '" << ns::kFDaily << "'!" << std::endl;
      return EXIT_FAILURE;
    }

    // "eop.txt" 書き込み
    if (!ns::write_data(data)) {
      std::cout << "[ERROR] Could not write '" << ns::kFEop << "'!" << std::endl;
      return EXIT_FAILURE;
    }
  } catch (...) {
      std::cerr << "EXCEPTION!" << std::endl;
      return EXIT_FAILURE;
  }

  return EXIT_SUCCESS;
}

• Gist - C++ source code to make an EOP(Earth Orientation Parameters) file.

3. C++ ソースコードのコンパイル

$ g++ -std=c++17 -Wall -O2 --pedantic-errors -o make_eop make_eop.cpp

4. 動作確認

$ ./make_eop

5. 結果確認

同じディレクトリ内に eop.txt というファイルが生成されるはずなので、確認してみる。

File: eop.txt

1973-01-02 41684.00 F  0.143000 0.000000  0.137000 0.000000 F  0.8075000 0.0000000  0.0000 0.1916 F   -18.637    0.000    -3.667    0.000
1973-01-03 41685.00 F  0.141000 0.000000  0.134000 0.000000 F  0.8044000 0.0000000  3.5563 0.1916 F   -18.636    0.000    -3.571    0.000
1973-01-04 41686.00 F  0.139000 0.000000  0.131000 0.000000 F  0.8012000 0.0000000  2.6599 0.1916 F   -18.669    0.000    -3.621    0.000
1973-01-05 41687.00 F  0.137000 0.000000  0.128000 0.000000 F  0.7981000 0.0000000  3.0344 0.1916 F   -18.751    0.000    -3.769    0.000
1973-01-06 41688.00 F  0.136000 0.000000  0.126000 0.000000 F  0.7949000 0.0000000  3.1276 0.1916 F   -18.868    0.000    -3.868    0.000
1973-01-07 41689.00 F  0.134000 0.000000  0.123000 0.000000 F  0.7918000 0.0000000  3.3271 0.1916 F   -18.959    0.000    -3.842    0.000
1973-01-08 41690.00 F  0.132000 0.000000  0.122000 0.000000 F  0.7887000 0.0000000  3.4940 0.1916 F   -18.976    0.000    -3.750    0.000
1973-01-09 41691.00 F  0.130000 0.000000  0.120000 0.000000 F  0.7855000 0.0000000  3.6275 0.1916 F   -18.934    0.000    -3.698    0.000
1973-01-10 41692.00 F  0.128000 0.000000  0.119000 0.000000 F  0.7824000 0.0000000  3.6924 0.2189 F   -18.889    0.000    -3.732    0.000
1973-01-11 41693.00 F  0.126000 0.000000  0.117000 0.000000 F  0.7792000 0.0000000  3.6729 0.2373 F   -18.861    0.000    -3.846    0.000

      :
===< 中略 >===
      :

2021-06-20 59385.00 P  0.189556 0.007482  0.424129 0.009675 P -0.1529675 0.0078428                                                       
2021-06-21 59386.00 P  0.190942 0.007534  0.423400 0.009759 P -0.1526717 0.0079197                                                       
2021-06-22 59387.00 P  0.192318 0.007586  0.422648 0.009844 P -0.1522055 0.0079964                                                       
2021-06-23 59388.00 P  0.193681 0.007637  0.421874 0.009928 P -0.1514893 0.0080728                                                       
2021-06-24 59389.00 P  0.195033 0.007688  0.421078 0.010012 P -0.1505174 0.0081489                                                       
2021-06-25 59390.00 P  0.196372 0.007739  0.420260 0.010096 P -0.1493570 0.0082248                                                       
2021-06-26 59391.00 P  0.197698 0.007789  0.419421 0.010179 P -0.1481211 0.0083005                                                       
2021-06-27 59392.00 P  0.199011 0.007840  0.418560 0.010262 P -0.1469268 0.0083759                                                       
2021-06-28 59393.00 P  0.200310 0.007890  0.417678 0.010345 P -0.1458554 0.0084511                                                       
2021-06-29 59394.00 P  0.201594 0.007940  0.416775 0.010428 P -0.1449344 0.0085261                                                       

---

以上。