code:
#include "HamlibTransceiver.hpp"

#include <cstring>
#include <cmath>

#include <QByteArray>
#include <QString>
#include <QStandardPaths>
#include <QFile>
#include <QJsonDocument>
#include <QJsonObject>
#include <QJsonValue>
#include <QDebug>

#include "moc_HamlibTransceiver.cpp"

namespace
{
  // Unfortunately bandwidth is conflated  with mode, this is probably
  // because Icom do  the same. So we have to  care about bandwidth if
  // we want  to set  mode otherwise we  will end up  setting unwanted
  // bandwidths every time we change mode.  The best we can do via the
  // Hamlib API is to request the  normal option for the mode and hope
  // that an appropriate filter is selected.  Also ensure that mode is
  // only set is absolutely necessary.  On Icoms (and probably others)
  // the filter is  selected by number without checking  the actual BW
  // so unless the  "normal" defaults are set on the  rig we won't get
  // desirable results.
  //
  // As  an ultimate  workaround make  sure  the user  always has  the
  // option to skip mode setting altogether.

  // reroute Hamlib diagnostic messages to Qt
  int debug_callback (enum rig_debug_level_e level, rig_ptr_t /* arg */, char const * format, va_list ap)
  {
    QString message;
    static char constexpr fmt[] = "Hamlib: %s";
    message = message.vsprintf (format, ap).trimmed ();

    switch (level)
      {
      case RIG_DEBUG_BUG:
        qFatal (fmt, message.toLocal8Bit ().data ());
        break;

      case RIG_DEBUG_ERR:
        qCritical (fmt, message.toLocal8Bit ().data ());
        break;

      case RIG_DEBUG_WARN:
        qWarning (fmt, message.toLocal8Bit ().data ());
        break;

      default:
        qDebug (fmt, message.toLocal8Bit ().data ());
        break;
      }

    return 0;
  }

  // callback function that receives transceiver capabilities from the
  // hamlib libraries
  int register_callback (rig_caps const * caps, void * callback_data)
  {
    TransceiverFactory::Transceivers * rigs = reinterpret_cast<TransceiverFactory::Transceivers *> (callback_data);

    QString key;
    if (RIG_MODEL_DUMMY == caps->rig_model)
      {
        key = TransceiverFactory::basic_transceiver_name_;
      }
    else
      {
        key = QString::fromLatin1 (caps->mfg_name).trimmed ()
          + ' '+ QString::fromLatin1 (caps->model_name).trimmed ()
          // + ' '+ QString::fromLatin1 (caps->version).trimmed ()
          // + " (" + QString::fromLatin1 (rig_strstatus (caps->status)).trimmed () + ')'
          ;
      }

    auto port_type = TransceiverFactory::Capabilities::none;
    switch (caps->port_type)
      {
      case RIG_PORT_SERIAL:
        port_type = TransceiverFactory::Capabilities::serial;
        break;

      case RIG_PORT_NETWORK:
        port_type = TransceiverFactory::Capabilities::network;
        break;

      case RIG_PORT_USB:
        port_type = TransceiverFactory::Capabilities::usb;
        break;

      default: break;
      }
    (*rigs)[key] = TransceiverFactory::Capabilities (caps->rig_model
                                                     , port_type
                                                     , RIG_MODEL_DUMMY != caps->rig_model
                                                     && (RIG_PTT_RIG == caps->ptt_type
                                                         || RIG_PTT_RIG_MICDATA == caps->ptt_type)
                                                     , RIG_PTT_RIG_MICDATA == caps->ptt_type);

    return 1;			// keep them coming
  }

  int unregister_callback (rig_caps const * caps, void *)
  {
    rig_unregister (caps->rig_model);
    return 1;			// keep them coming
  }

  // int frequency_change_callback (RIG * /* rig */, vfo_t vfo, freq_t f, rig_ptr_t arg)
  // {
  //   (void)vfo;			// unused in release build

  //   Q_ASSERT (vfo == RIG_VFO_CURR); // G4WJS: at the time of writing only current VFO is signalled by hamlib

  //   HamlibTransceiver * transceiver (reinterpret_cast<HamlibTransceiver *> (arg));
  //   Q_EMIT transceiver->frequency_change (f, Transceiver::A);
  //   return RIG_OK;
  // }

  class hamlib_tx_vfo_fixup final
  {
  public:
    hamlib_tx_vfo_fixup (RIG * rig, vfo_t tx_vfo)
      : rig_ {rig}
    {
      original_vfo_ = rig_->state.tx_vfo;
      rig_->state.tx_vfo = tx_vfo;
    }

    ~hamlib_tx_vfo_fixup ()
    {
      rig_->state.tx_vfo = original_vfo_;
    }

  private:
    RIG * rig_;
    vfo_t original_vfo_;
  };
}

freq_t HamlibTransceiver::dummy_frequency_;
rmode_t HamlibTransceiver::dummy_mode_ {RIG_MODE_NONE};

void HamlibTransceiver::register_transceivers (TransceiverFactory::Transceivers * registry)
{
  rig_set_debug_callback (debug_callback, nullptr);

#if WSJT_HAMLIB_TRACE
#if WSJT_HAMLIB_VERBOSE_TRACE
  rig_set_debug (RIG_DEBUG_TRACE);
#else
  rig_set_debug (RIG_DEBUG_VERBOSE);
#endif
#elif defined (NDEBUG)
  rig_set_debug (RIG_DEBUG_ERR);
#else
  rig_set_debug (RIG_DEBUG_WARN);
#endif

  rig_load_all_backends ();
  rig_list_foreach (register_callback, registry);
}

void HamlibTransceiver::unregister_transceivers ()
{
  rig_list_foreach (unregister_callback, nullptr);
}

void HamlibTransceiver::RIGDeleter::cleanup (RIG * rig)
{
  if (rig)
    {
      // rig->state.obj = 0;
      rig_cleanup (rig);
    }
}

HamlibTransceiver::HamlibTransceiver (TransceiverFactory::PTTMethod ptt_type, QString const& ptt_port,
                                      QObject * parent)
  : PollingTransceiver {0, parent}
  , rig_ {rig_init (RIG_MODEL_DUMMY)}
  , back_ptt_port_ {false}
  , one_VFO_ {false}
  , is_dummy_ {true}
  , reversed_ {false}
  , freq_query_works_ {true}
  , mode_query_works_ {true}
  , split_query_works_ {true}
  , tickle_hamlib_ {false}
  , get_vfo_works_ {true}
  , set_vfo_works_ {true}
{
  if (!rig_)
    {
      throw error {tr ("Hamlib initialisation error")};
    }

  switch (ptt_type)
    {
    case TransceiverFactory::PTT_method_VOX:
      set_conf ("ptt_type", "None");
      break;

    case TransceiverFactory::PTT_method_CAT:
      // Use the default PTT_TYPE for the rig (defined in the Hamlib
      // rig back-end capabilities).
      break;

    case TransceiverFactory::PTT_method_DTR:
    case TransceiverFactory::PTT_method_RTS:
      if (!ptt_port.isEmpty ())
        {
#if defined (WIN32)
          set_conf ("ptt_pathname", ("\\\\.\\" + ptt_port).toLatin1 ().data ());
#else
          set_conf ("ptt_pathname", ptt_port.toLatin1 ().data ());
#endif
        }

      if (TransceiverFactory::PTT_method_DTR == ptt_type)
        {
          set_conf ("ptt_type", "DTR");
        }
      else
        {
          set_conf ("ptt_type", "RTS");
        }
    }
}

HamlibTransceiver::HamlibTransceiver (int model_number, TransceiverFactory::ParameterPack const& params,
                                      QObject * parent)
  : PollingTransceiver {params.poll_interval, parent}
  , rig_ {rig_init (model_number)}
  , back_ptt_port_ {TransceiverFactory::TX_audio_source_rear == params.audio_source}
  , one_VFO_ {false}
  , is_dummy_ {RIG_MODEL_DUMMY == model_number}
  , reversed_ {false}
  , freq_query_works_ {rig_ && rig_->caps->get_freq}
  , mode_query_works_ {rig_ && rig_->caps->get_mode}
  , split_query_works_ {rig_ && rig_->caps->get_split_vfo}
  , tickle_hamlib_ {false}
  , get_vfo_works_ {true}
  , set_vfo_works_ {true}
{
  if (!rig_)
    {
      throw error {tr ("Hamlib initialisation error")};
    }

  // rig_->state.obj = this;

  if (!is_dummy_)
    {
      //
      // user defined Hamlib settings
      //
      auto settings_file_name = QStandardPaths::locate (
#if QT_VERSION >= 0x050500
                                                        QStandardPaths::AppConfigLocation
#else
                                                        QStandardPaths::ConfigLocation
#endif
                                                        , "hamlib_settings.json");
      if (!settings_file_name.isEmpty ())
        {
          QFile settings_file {settings_file_name};
          qDebug () << "Using Hamlib settings file:" << settings_file_name;
          if (settings_file.open (QFile::ReadOnly))
            {
              QJsonParseError status;
              auto settings_doc = QJsonDocument::fromJson (settings_file.readAll (), &status);
              if (status.error)
                {
                  throw error {tr ("Hamlib settings file error: %1 at character offset %2")
                      .arg (status.errorString ()).arg (status.offset)};
                }
              qDebug () << "Hamlib settings JSON:" << settings_doc.toJson ();
              if (!settings_doc.isObject ())
                {
                  throw error {tr ("Hamlib settings file error: top level must be a JSON object")};
                }
              auto const& settings = settings_doc.object ();

              //
              // configuration settings
              //
              auto const& config = settings["config"];
              if (!config.isUndefined ())
                {
                  if (!config.isObject ())
                    {
                      throw error {tr ("Hamlib settings file error: config must be a JSON object")};
                    }
                  auto const& config_list = config.toObject ();
                  for (auto item = config_list.constBegin (); item != config_list.constEnd (); ++item)
                    {
                      set_conf (item.key ().toLocal8Bit ().constData ()
                                , (*item).toVariant ().toString ().toLocal8Bit ().constData ());
                    }
                }
            }
        }

      switch (rig_->caps->port_type)
        {
        case RIG_PORT_SERIAL:
          if (!params.serial_port.isEmpty ())
            {
              set_conf ("rig_pathname", params.serial_port.toLatin1 ().data ());
            }
          set_conf ("serial_speed", QByteArray::number (params.baud).data ());
          if (params.data_bits != TransceiverFactory::default_data_bits)
            {
              set_conf ("data_bits", TransceiverFactory::seven_data_bits == params.data_bits ? "7" : "8");
            }
          if (params.stop_bits != TransceiverFactory::default_stop_bits)
            {
              set_conf ("stop_bits", TransceiverFactory::one_stop_bit == params.stop_bits ? "1" : "2");
            }

          switch (params.handshake)
            {
            case TransceiverFactory::handshake_none: set_conf ("serial_handshake", "None"); break;
            case TransceiverFactory::handshake_XonXoff: set_conf ("serial_handshake", "XONXOFF"); break;
            case TransceiverFactory::handshake_hardware: set_conf ("serial_handshake", "Hardware"); break;
            default: break;
            }

          if (params.force_dtr)
            {
              set_conf ("dtr_state", params.dtr_high ? "ON" : "OFF");
            }
          if (params.force_rts)
            {
              if (TransceiverFactory::handshake_hardware != params.handshake)
                {
                  set_conf ("rts_state", params.rts_high ? "ON" : "OFF");
                }
            }
          break;

        case RIG_PORT_NETWORK:
          if (!params.network_port.isEmpty ())
            {
              set_conf ("rig_pathname", params.network_port.toLatin1 ().data ());
            }
          break;

        case RIG_PORT_USB:
          if (!params.usb_port.isEmpty ())
            {
              set_conf ("rig_pathname", params.usb_port.toLatin1 ().data ());
            }
          break;

        default:
          throw error {tr ("Unsupported CAT type")};
          break;
        }
    }

  switch (params.ptt_type)
    {
    case TransceiverFactory::PTT_method_VOX:
      set_conf ("ptt_type", "None");
      break;

    case TransceiverFactory::PTT_method_CAT:
      // Use the default PTT_TYPE for the rig (defined in the Hamlib
      // rig back-end capabilities).
      break;

    case TransceiverFactory::PTT_method_DTR:
    case TransceiverFactory::PTT_method_RTS:
      if (!params.ptt_port.isEmpty ()
          && params.ptt_port != "None"
          && (is_dummy_ || params.ptt_port != params.serial_port))
        {
#if defined (WIN32)
          set_conf ("ptt_pathname", ("\\\\.\\" + params.ptt_port).toLatin1 ().data ());
#else
          set_conf ("ptt_pathname", params.ptt_port.toLatin1 ().data ());
#endif
        }

      if (TransceiverFactory::PTT_method_DTR == params.ptt_type)
        {
          set_conf ("ptt_type", "DTR");
        }
      else
        {
          set_conf ("ptt_type", "RTS");
        }
    }

  // Make Icom CAT split commands less glitchy
  set_conf ("no_xchg", "1");

  // would be nice to get events but not supported on Windows and also not on a lot of rigs
  // rig_set_freq_callback (rig_.data (), &frequency_change_callback, this);
}

void HamlibTransceiver::error_check (int ret_code, QString const& doing) const
{
  if (RIG_OK != ret_code)
    {
      TRACE_CAT_POLL ("HamlibTransceiver", "error:" << rigerror (ret_code));
      throw error {tr ("Hamlib error: %1 while %2").arg (rigerror (ret_code)).arg (doing)};
    }
}

int HamlibTransceiver::do_start ()
{
  TRACE_CAT ("HamlibTransceiver",
             QString::fromLatin1 (rig_->caps->mfg_name).trimmed ()
             << QString::fromLatin1 (rig_->caps->model_name).trimmed ());

  error_check (rig_open (rig_.data ()), tr ("opening connection to rig"));

  // reset dynamic state
  one_VFO_ = false;
  reversed_ = false;
  freq_query_works_ = rig_->caps->get_freq;
  mode_query_works_ = rig_->caps->get_mode;
  split_query_works_ = rig_->caps->get_split_vfo;
  tickle_hamlib_ = false;
  get_vfo_works_ = true;
  set_vfo_works_ = true;

  // the Net rigctl back end promises all functions work but we must
  // test get_vfo as it determines our strategy for Icom rigs
  vfo_t vfo;
  int rc = rig_get_vfo (rig_.data (), &vfo);
  if (-RIG_ENAVAIL == rc || -RIG_ENIMPL == rc)
    {
      get_vfo_works_ = false;
      // determine if the rig uses single VFO addressing i.e. A/B and
      // no get_vfo function
      if (rig_->state.vfo_list & RIG_VFO_B)
        {
          one_VFO_ = true;
        }
    }
  else
    {
      error_check (rc, "testing getting current VFO");
    }

  if ((WSJT_RIG_NONE_CAN_SPLIT || !is_dummy_)
      && rig_->caps->set_split_vfo) // if split is possible do some extra setup
    {
      freq_t f1;
      freq_t f2;
      rmode_t m {RIG_MODE_USB};
      rmode_t mb;
      pbwidth_t w {RIG_PASSBAND_NORMAL};
      pbwidth_t wb;
      if (freq_query_works_
          && (!get_vfo_works_ || !rig_->caps->get_vfo))
        {
          // Icom have deficient CAT protocol with no way of reading which
          // VFO is selected or if SPLIT is selected so we have to simply
          // assume it is as when we started by setting at open time right
          // here. We also gather/set other initial state.
          error_check (rig_get_freq (rig_.data (), RIG_VFO_CURR, &f1), tr ("getting current frequency"));
          f1 = std::round (f1);
          TRACE_CAT ("HamlibTransceiver", "current frequency =" << f1);

          error_check (rig_get_mode (rig_.data (), RIG_VFO_CURR, &m, &w), tr ("getting current mode"));
          TRACE_CAT ("HamlibTransceiver", "current mode =" << rig_strrmode (m) << "bw =" << w);

          if (!rig_->caps->set_vfo)
            {
              TRACE_CAT ("HamlibTransceiver", "rig_vfo_op TOGGLE");
              rc = rig_vfo_op (rig_.data (), RIG_VFO_CURR, RIG_OP_TOGGLE);
            }
          else
            {
              TRACE_CAT ("HamlibTransceiver", "rig_set_vfo to other VFO");
              rc = rig_set_vfo (rig_.data (), rig_->state.vfo_list & RIG_VFO_B ? RIG_VFO_B : RIG_VFO_SUB);
              if (-RIG_ENAVAIL == rc || -RIG_ENIMPL == rc)
                {
                  // if we are talking to netrigctl then toggle VFO op
                  // may still work
                  TRACE_CAT ("HamlibTransceiver", "rig_vfo_op TOGGLE");
                  rc = rig_vfo_op (rig_.data (), RIG_VFO_CURR, RIG_OP_TOGGLE);
                }
            }
          if (-RIG_ENAVAIL == rc || -RIG_ENIMPL == rc)
            {
              // we are probably dealing with rigctld so we do not
              // have completely accurate rig capabilities
              set_vfo_works_ = false;
              one_VFO_ = false; // we do not need single VFO addressing
            }
          else
            {
              error_check (rc, tr ("exchanging VFOs"));
            }

          if (set_vfo_works_)
            {
              // without the above we cannot proceed but we know we
              // are on VFO A and that will not change so there's no
              // need to execute this block
              error_check (rig_get_freq (rig_.data (), RIG_VFO_CURR, &f2), tr ("getting other VFO frequency"));
              f2 = std::round (f2);
              TRACE_CAT ("HamlibTransceiver", "rig_get_freq other frequency =" << f2);

              error_check (rig_get_mode (rig_.data (), RIG_VFO_CURR, &mb, &wb), tr ("getting other VFO mode"));
              TRACE_CAT ("HamlibTransceiver", "rig_get_mode other mode =" << rig_strrmode (mb) << "bw =" << wb);

              update_other_frequency (f2);

              if (!rig_->caps->set_vfo)
                {
                  TRACE_CAT ("HamlibTransceiver", "rig_vfo_op TOGGLE");
                  error_check (rig_vfo_op (rig_.data (), RIG_VFO_CURR, RIG_OP_TOGGLE), tr ("exchanging VFOs"));
                }
              else
                {
                  TRACE_CAT ("HamlibTransceiver", "rig_set_vfo A/MAIN");
                  error_check (rig_set_vfo (rig_.data (), rig_->state.vfo_list & RIG_VFO_A ? RIG_VFO_A : RIG_VFO_MAIN), tr ("setting current VFO"));
                }

              if (f1 != f2 || m != mb || w != wb)	// we must have started with MAIN/A
                {
                  update_rx_frequency (f1);
                }
              else
                {
                  error_check (rig_get_freq (rig_.data (), RIG_VFO_CURR, &f1), tr ("getting frequency"));
                  f1 = std::round (f1);
                  TRACE_CAT ("HamlibTransceiver", "rig_get_freq frequency =" << f1);

                  error_check (rig_get_mode (rig_.data (), RIG_VFO_CURR, &m, &w), tr ("getting mode"));
                  TRACE_CAT ("HamlibTransceiver", "rig_get_mode mode =" << rig_strrmode (m) << "bw =" << w);

                  update_rx_frequency (f1);
                }
            }

          // TRACE_CAT ("HamlibTransceiver", "rig_set_split_vfo split off");
          // error_check (rig_set_split_vfo (rig_.data (), RIG_VFO_CURR, RIG_SPLIT_OFF, RIG_VFO_CURR), tr ("setting split off"));
          // update_split (false);
        }
      else
        {
          vfo_t v {RIG_VFO_A};  // assume RX always on VFO A/MAIN

          if (get_vfo_works_ && rig_->caps->get_vfo)
            {
              error_check (rig_get_vfo (rig_.data (), &v), tr ("getting current VFO")); // has side effect of establishing current VFO inside hamlib
              TRACE_CAT ("HamlibTransceiver", "rig_get_vfo current VFO = " << rig_strvfo (v));
            }

          reversed_ = RIG_VFO_B == v;

          if (mode_query_works_ && !(rig_->caps->targetable_vfo & (RIG_TARGETABLE_MODE | RIG_TARGETABLE_PURE)))
            {
              if (RIG_OK == rig_get_mode (rig_.data (), RIG_VFO_CURR, &m, &w))
                {
                  TRACE_CAT ("HamlibTransceiver", "rig_get_mode current mode =" << rig_strrmode (m) << "bw =" << w);
                }
              else
                {
                  mode_query_works_ = false;
                  // Some rigs (HDSDR) don't have a working way of
                  // reporting MODE so we give up on mode queries -
                  // sets will still cause an error
                  TRACE_CAT ("HamlibTransceiver", "rig_get_mode can't do on this rig");
                }
            }
        }
      update_mode (map_mode (m));
    }

  tickle_hamlib_ = true;

  if (is_dummy_ && dummy_frequency_)
    {
      // return to where last dummy instance was
      // TODO: this is going to break down if multiple dummy rigs are used
      rig_set_freq (rig_.data (), RIG_VFO_CURR, dummy_frequency_);
      update_rx_frequency (dummy_frequency_);
      if (RIG_MODE_NONE != dummy_mode_)
        {
          rig_set_mode (rig_.data (), RIG_VFO_CURR, dummy_mode_, RIG_PASSBAND_NOCHANGE);
          update_mode (map_mode (dummy_mode_));
        }
    }

  int resolution {0};
  if (freq_query_works_)
    {
      freq_t current_frequency;
      error_check (rig_get_freq (rig_.data (), RIG_VFO_CURR, &current_frequency), tr ("getting current VFO frequency"));
      current_frequency = std::round (current_frequency);
      Frequency f = current_frequency;
      if (f && !(f % 10))
        {
          auto test_frequency = f - f % 100 + 55;
          error_check (rig_set_freq (rig_.data (), RIG_VFO_CURR, test_frequency), tr ("setting frequency"));
          freq_t new_frequency;
          error_check (rig_get_freq (rig_.data (), RIG_VFO_CURR, &new_frequency), tr ("getting current VFO frequency"));
          new_frequency = std::round (new_frequency);
          switch (static_cast<Radio::FrequencyDelta> (new_frequency - test_frequency))
            {
            case -5: resolution = -1; break;  // 10Hz truncated
            case 5: resolution = 1; break;    // 10Hz rounded
            case -15: resolution = -2; break; // 20Hz truncated
            case -55: resolution = -3; break; // 100Hz truncated
            case 45: resolution = 3; break;   // 100Hz rounded
            }
          if (1 == resolution)      // may be 20Hz rounded
            {
              test_frequency = f - f % 100 + 51;
              error_check (rig_set_freq (rig_.data (), RIG_VFO_CURR, test_frequency), tr ("setting frequency"));
              error_check (rig_get_freq (rig_.data (), RIG_VFO_CURR, &new_frequency), tr ("getting current VFO frequency"));
              if (9 == static_cast<Radio::FrequencyDelta> (new_frequency - test_frequency))
                {
                  resolution = 2;   // 20Hz rounded
                }
            }
          error_check (rig_set_freq (rig_.data (), RIG_VFO_CURR, current_frequency), tr ("setting frequency"));
        }
    }
  else
    {
      resolution = -1;          // best guess
    }

  poll ();

  TRACE_CAT ("HamlibTransceiver", "exit" << state () << "reversed =" << reversed_ << "resolution = " << resolution);
  return resolution;
}

void HamlibTransceiver::do_stop ()
{
  if (is_dummy_)
    {
      rig_get_freq (rig_.data (), RIG_VFO_CURR, &dummy_frequency_);
      dummy_frequency_ = std::round (dummy_frequency_);
      if (mode_query_works_)
        {
          pbwidth_t width;
          rig_get_mode (rig_.data (), RIG_VFO_CURR, &dummy_mode_, &width);
        }
    }
  if (rig_)
    {
      rig_close (rig_.data ());
    }

  TRACE_CAT ("HamlibTransceiver", "state:" << state () << "reversed =" << reversed_);
}

auto HamlibTransceiver::get_vfos (bool for_split) const -> std::tuple<vfo_t, vfo_t>
{
  if (get_vfo_works_ && rig_->caps->get_vfo)
    {
      vfo_t v;
      error_check (rig_get_vfo (rig_.data (), &v), tr ("getting current VFO")); // has side effect of establishing current VFO inside hamlib
      TRACE_CAT ("HamlibTransceiver", "rig_get_vfo VFO = " << rig_strvfo (v));

      reversed_ = RIG_VFO_B == v;
    }
  else if (!for_split && set_vfo_works_ && rig_->caps->set_vfo && rig_->caps->set_split_vfo)
    {
      // use VFO A/MAIN for main frequency and B/SUB for Tx
      // frequency if split since these type of radios can only
      // support this way around

      TRACE_CAT ("HamlibTransceiver", "rig_set_vfo VFO = A/MAIN");
      error_check (rig_set_vfo (rig_.data (), rig_->state.vfo_list & RIG_VFO_A ? RIG_VFO_A : RIG_VFO_MAIN), tr ("setting current VFO"));
    }
  // else only toggle available but VFOs should be substitutable

  auto rx_vfo = rig_->state.vfo_list & RIG_VFO_A ? RIG_VFO_A : RIG_VFO_MAIN;
  auto tx_vfo = (WSJT_RIG_NONE_CAN_SPLIT || !is_dummy_) && for_split
    ? (rig_->state.vfo_list & RIG_VFO_B ? RIG_VFO_B : RIG_VFO_SUB)
    : rx_vfo;
  if (reversed_)
    {
      TRACE_CAT ("HamlibTransceiver", "reversing VFOs");
      std::swap (rx_vfo, tx_vfo);
    }

  TRACE_CAT ("HamlibTransceiver", "RX VFO = " << rig_strvfo (rx_vfo) << " TX VFO = " << rig_strvfo (tx_vfo));
  return std::make_tuple (rx_vfo, tx_vfo);
}

void HamlibTransceiver::do_frequency (Frequency f, MODE m, bool no_ignore)
{
  TRACE_CAT ("HamlibTransceiver", f << "mode:" << m << "reversed:" << reversed_);

  // only change when receiving or simplex or direct VFO addressing
  // unavailable or forced
  if (!state ().ptt () || !state ().split () || !one_VFO_ || no_ignore)
    {
      // for the 1st time as a band change may cause a recalled mode to be
      // set
      error_check (rig_set_freq (rig_.data (), RIG_VFO_CURR, f), tr ("setting frequency"));
      update_rx_frequency (f);

      if (mode_query_works_ && UNK != m)
        {
          rmode_t current_mode;
          pbwidth_t current_width;
          auto new_mode = map_mode (m);
          error_check (rig_get_mode (rig_.data (), RIG_VFO_CURR, &current_mode, &current_width), tr ("getting current VFO mode"));
          TRACE_CAT ("HamlibTransceiver", "rig_get_mode mode = " << rig_strrmode (current_mode) << "bw =" << current_width);

          if (new_mode != current_mode)
            {
              TRACE_CAT ("HamlibTransceiver", "rig_set_mode mode = " << rig_strrmode (new_mode));
              error_check (rig_set_mode (rig_.data (), RIG_VFO_CURR, new_mode, RIG_PASSBAND_NOCHANGE), tr ("setting current VFO mode"));

              // for the 2nd time because a mode change may have caused a
              // frequency change
              error_check (rig_set_freq (rig_.data (), RIG_VFO_CURR, f), tr ("setting frequency"));

              // for the second time because some rigs change mode according
              // to frequency such as the TS-2000 auto mode setting
              TRACE_CAT ("HamlibTransceiver", "rig_set_mode mode = " << rig_strrmode (new_mode));
              error_check (rig_set_mode (rig_.data (), RIG_VFO_CURR, new_mode, RIG_PASSBAND_NOCHANGE), tr ("setting current VFO mode"));
            }
          update_mode (m);
        }
    }
}

void HamlibTransceiver::do_tx_frequency (Frequency tx, MODE mode, bool no_ignore)
{
  TRACE_CAT ("HamlibTransceiver", tx << "reversed:" << reversed_);

  if (WSJT_RIG_NONE_CAN_SPLIT || !is_dummy_) // split is meaningless if you can't see it
    {
      auto split = tx ? RIG_SPLIT_ON : RIG_SPLIT_OFF;
      auto vfos = get_vfos (tx);
      // auto rx_vfo = std::get<0> (vfos); // or use RIG_VFO_CURR
      auto tx_vfo = std::get<1> (vfos);

      if (tx)
        {
          // Doing set split for the 1st of two times, this one
          // ensures that the internal Hamlib state is correct
          // otherwise rig_set_split_freq() will target the wrong VFO
          // on some rigs

          if (tickle_hamlib_)
            {
              // This potentially causes issues with the Elecraft K3
              // which will block setting split mode when it deems
              // cross mode split operation not possible. There's not
              // much we can do since the Hamlib Library needs this
              // call at least once to establish the Tx VFO. Best we
              // can do is only do this once per session.
              TRACE_CAT ("HamlibTransceiver", "rig_set_split_vfo split =" << split);
              auto rc = rig_set_split_vfo (rig_.data (), RIG_VFO_CURR, split, tx_vfo);
              if (tx || (-RIG_ENAVAIL != rc && -RIG_ENIMPL != rc))
                {
                  // On rigs that can't have split controlled only throw an
                  // exception when an error other than command not accepted
                  // is returned when trying to leave split mode. This allows
                  // fake split mode and non-split mode to work without error
                  // on such rigs without having to know anything about the
                  // specific rig.
                  error_check (rc, tr ("setting/unsetting split mode"));
                }
              tickle_hamlib_ = false;
              update_split (tx);
            }

          // just change current when transmitting with single VFO
          // addressing
          if (state ().ptt () && one_VFO_)
            {
              TRACE_CAT ("HamlibTransceiver", "rig_set_split_vfo split =" << split);
              error_check (rig_set_split_vfo (rig_.data (), RIG_VFO_CURR, split, tx_vfo), tr ("setting split mode"));

              error_check (rig_set_freq (rig_.data (), RIG_VFO_CURR, tx), tr ("setting frequency"));

              if (UNK != mode && mode_query_works_)
                {
                  rmode_t current_mode;
                  pbwidth_t current_width;
                  auto new_mode = map_mode (mode);
                  error_check (rig_get_mode (rig_.data (), RIG_VFO_CURR, &current_mode, &current_width), tr ("getting current VFO mode"));
                  TRACE_CAT ("HamlibTransceiver", "rig_get_mode mode = " << rig_strrmode (current_mode) << "bw =" << current_width);

                  if (new_mode != current_mode)
                    {
                      TRACE_CAT ("HamlibTransceiver", "rig_set_mode mode = " << rig_strrmode (new_mode));
                      error_check (rig_set_mode (rig_.data (), RIG_VFO_CURR, new_mode, RIG_PASSBAND_NOCHANGE), tr ("setting current VFO mode"));
                    }
                }
              update_other_frequency (tx);
            }
          else if (!one_VFO_ || no_ignore)   // if not single VFO addressing and not forced
            {
              hamlib_tx_vfo_fixup fixup (rig_.data (), tx_vfo);
              if (UNK != mode)
                {
                  auto new_mode = map_mode (mode);
                  TRACE_CAT ("HamlibTransceiver", "rig_set_split_freq_mode freq = " << tx
                             << " mode = " << rig_strrmode (new_mode));
                  error_check (rig_set_split_freq_mode (rig_.data (), RIG_VFO_CURR, tx, new_mode, RIG_PASSBAND_NOCHANGE), tr ("setting split TX frequency and mode"));
                }
              else
                {
                  TRACE_CAT ("HamlibTransceiver", "rig_set_split_freq freq = " << tx);
                  error_check (rig_set_split_freq (rig_.data (), RIG_VFO_CURR, tx), tr ("setting split TX frequency"));
                }
              // Enable split last since some rigs (Kenwood for one) come out
              // of split when you switch RX VFO (to set split mode above for
              // example). Also the Elecraft K3 will refuse to go to split
              // with certain VFO A/B mode combinations.
              TRACE_CAT ("HamlibTransceiver", "rig_set_split_vfo split =" << split);
              error_check (rig_set_split_vfo (rig_.data (), RIG_VFO_CURR, split, tx_vfo), tr ("setting split mode"));
              update_other_frequency (tx);
              update_split (tx);
            }
        }
      else
        {
          // Disable split
          TRACE_CAT ("HamlibTransceiver", "rig_set_split_vfo split =" << split);
          auto rc = rig_set_split_vfo (rig_.data (), RIG_VFO_CURR, split, tx_vfo);
          if (tx || (-RIG_ENAVAIL != rc && -RIG_ENIMPL != rc))
            {
              // On rigs that can't have split controlled only throw an
              // exception when an error other than command not accepted
              // is returned when trying to leave split mode. This allows
              // fake split mode and non-split mode to work without error
              // on such rigs without having to know anything about the
              // specific rig.
              error_check (rc, tr ("setting/unsetting split mode"));
            }
          update_other_frequency (tx);
          update_split (tx);
        }
    }
}

void HamlibTransceiver::do_mode (MODE mode)
{
  TRACE_CAT ("HamlibTransceiver", mode);

  auto vfos = get_vfos (state ().split ());
  // auto rx_vfo = std::get<0> (vfos);
  auto tx_vfo = std::get<1> (vfos);

  rmode_t current_mode;
  pbwidth_t current_width;
  auto new_mode = map_mode (mode);

  // only change when receiving or simplex if direct VFO addressing unavailable
  if (!(state ().ptt () && state ().split () && one_VFO_))
    {
      error_check (rig_get_mode (rig_.data (), RIG_VFO_CURR, &current_mode, &current_width), tr ("getting current VFO mode"));
      TRACE_CAT ("HamlibTransceiver", "rig_get_mode mode = " << rig_strrmode (current_mode) << "bw =" << current_width);

      if (new_mode != current_mode)
        {
          TRACE_CAT ("HamlibTransceiver", "rig_set_mode mode = " << rig_strrmode (new_mode));
          error_check (rig_set_mode (rig_.data (), RIG_VFO_CURR, new_mode, RIG_PASSBAND_NOCHANGE), tr ("setting current VFO mode"));
        }
    }

  // just change current when transmitting split with one VFO mode
  if (state ().ptt () && state ().split () && one_VFO_)
    {
      error_check (rig_get_mode (rig_.data (), RIG_VFO_CURR, &current_mode, &current_width), tr ("getting current VFO mode"));
      TRACE_CAT ("HamlibTransceiver", "rig_get_mode mode = " << rig_strrmode (current_mode) << "bw =" << current_width);

      if (new_mode != current_mode)
        {
          TRACE_CAT ("HamlibTransceiver", "rig_set_mode mode = " << rig_strrmode (new_mode));
          error_check (rig_set_mode (rig_.data (), RIG_VFO_CURR, new_mode, RIG_PASSBAND_NOCHANGE), tr ("setting current VFO mode"));
        }
    }
  else if (state ().split () && !one_VFO_)
    {
      error_check (rig_get_split_mode (rig_.data (), RIG_VFO_CURR, &current_mode, &current_width), tr ("getting split TX VFO mode"));
      TRACE_CAT ("HamlibTransceiver", "rig_get_split_mode mode = " << rig_strrmode (current_mode) << "bw =" << current_width);

      if (new_mode != current_mode)
        {
          TRACE_CAT ("HamlibTransceiver", "rig_set_split_mode mode = " << rig_strrmode (new_mode));
          hamlib_tx_vfo_fixup fixup (rig_.data (), tx_vfo);
          error_check (rig_set_split_mode (rig_.data (), RIG_VFO_CURR, new_mode, RIG_PASSBAND_NOCHANGE), tr ("setting split TX VFO mode"));
        }
    }
  update_mode (mode);
}

void HamlibTransceiver::poll ()
{
#if !WSJT_TRACE_CAT_POLLS
#if defined (NDEBUG)
  rig_set_debug (RIG_DEBUG_ERR);
#else
  rig_set_debug (RIG_DEBUG_WARN);
#endif
#endif

  freq_t f;
  rmode_t m;
  pbwidth_t w;
  split_t s;

  if (get_vfo_works_ && rig_->caps->get_vfo)
    {
      vfo_t v;
      error_check (rig_get_vfo (rig_.data (), &v), tr ("getting current VFO")); // has side effect of establishing current VFO inside hamlib
      TRACE_CAT_POLL ("HamlibTransceiver", "VFO =" << rig_strvfo (v));
      reversed_ = RIG_VFO_B == v;
    }

  if ((WSJT_RIG_NONE_CAN_SPLIT || !is_dummy_)
      && rig_->caps->get_split_vfo && split_query_works_)
    {
      vfo_t v {RIG_VFO_NONE};		// so we can tell if it doesn't get updated :(
      auto rc = rig_get_split_vfo (rig_.data (), RIG_VFO_CURR, &s, &v);
      if (-RIG_OK == rc && RIG_SPLIT_ON == s)
        {
          TRACE_CAT_POLL ("HamlibTransceiver", "rig_get_split_vfo split = " << s << " VFO = " << rig_strvfo (v));
          update_split (true);
          // if (RIG_VFO_A == v)
          // 	{
          // 	  reversed_ = true;	// not sure if this helps us here
          // 	}
        }
      else if (-RIG_OK == rc)	// not split
        {
          TRACE_CAT_POLL ("HamlibTransceiver", "rig_get_split_vfo split = " << s << " VFO = " << rig_strvfo (v));
          update_split (false);
        }
      else
        {
          // Some rigs (Icom) don't have a way of reporting SPLIT
          // mode
          TRACE_CAT_POLL ("HamlibTransceiver", "rig_get_split_vfo can't do on this rig");
          // just report how we see it based on prior commands
          split_query_works_ = false;
        }
    }

  if (freq_query_works_)
    {
      // only read if possible and when receiving or simplex
      if (!state ().ptt () || !state ().split ())
        {
          error_check (rig_get_freq (rig_.data (), RIG_VFO_CURR, &f), tr ("getting current VFO frequency"));
          f = std::round (f);
          TRACE_CAT_POLL ("HamlibTransceiver", "rig_get_freq frequency =" << f);
          update_rx_frequency (f);
        }

      if ((WSJT_RIG_NONE_CAN_SPLIT || !is_dummy_)
          && state ().split ()
          && (rig_->caps->targetable_vfo & (RIG_TARGETABLE_FREQ | RIG_TARGETABLE_PURE))
          && !one_VFO_)
        {
          // only read "other" VFO if in split, this allows rigs like
          // FlexRadio to work in Kenwood TS-2000 mode despite them
          // not having a FB; command

          // we can only probe current VFO unless rig supports reading
          // the other one directly because we can't glitch the Rx
          error_check (rig_get_freq (rig_.data ()
                                     , reversed_
                                     ? (rig_->state.vfo_list & RIG_VFO_A ? RIG_VFO_A : RIG_VFO_MAIN)
                                     : (rig_->state.vfo_list & RIG_VFO_B ? RIG_VFO_B : RIG_VFO_SUB)
                                     , &f), tr ("getting other VFO frequency"));
          f = std::round (f);
          TRACE_CAT_POLL ("HamlibTransceiver", "rig_get_freq other VFO =" << f);
          update_other_frequency (f);
        }
    }

  // only read when receiving or simplex if direct VFO addressing unavailable
  if ((!state ().ptt () || !state ().split ())
      && mode_query_works_)
    {
      // We have to ignore errors here because Yaesu FTdx... rigs can
      // report the wrong mode when transmitting split with different
      // modes per VFO. This is unfortunate because that is exactly
      // what you need to do to get 4kHz Rx b.w and modulation into
      // the rig through the data socket or USB. I.e.  USB for Rx and
      // DATA-USB for Tx.
      auto rc = rig_get_mode (rig_.data (), RIG_VFO_CURR, &m, &w);
      if (RIG_OK == rc)
        {
          TRACE_CAT_POLL ("HamlibTransceiver", "rig_get_mode mode =" << rig_strrmode (m) << "bw =" << w);
          update_mode (map_mode (m));
        }
      else
        {
          TRACE_CAT_POLL ("HamlibTransceiver", "rig_get_mode mode failed with rc:" << rc << "ignoring");
        }
    }

  if (RIG_PTT_NONE != rig_->state.pttport.type.ptt && rig_->caps->get_ptt)
    {
      ptt_t p;
      auto rc = rig_get_ptt (rig_.data (), RIG_VFO_CURR, &p);
      if (-RIG_ENAVAIL != rc && -RIG_ENIMPL != rc) // may fail if
        // Net rig ctl and target doesn't
        // support command
        {
          error_check (rc, tr ("getting PTT state"));
          TRACE_CAT_POLL ("HamlibTransceiver", "rig_get_ptt PTT =" << p);
          update_PTT (!(RIG_PTT_OFF == p));
        }
    }

#if !WSJT_TRACE_CAT_POLLS
#if WSJT_HAMLIB_TRACE
#if WSJT_HAMLIB_VERBOSE_TRACE
  rig_set_debug (RIG_DEBUG_TRACE);
#else
  rig_set_debug (RIG_DEBUG_VERBOSE);
#endif
#elif defined (NDEBUG)
  rig_set_debug (RIG_DEBUG_ERR);
#else
  rig_set_debug (RIG_DEBUG_WARN);
#endif
#endif
}

void HamlibTransceiver::do_ptt (bool on)
{
  TRACE_CAT ("HamlibTransceiver", on << state () << "reversed =" << reversed_);
  if (on)
    {
      if (RIG_PTT_NONE != rig_->state.pttport.type.ptt)
        {
          TRACE_CAT ("HamlibTransceiver", "rig_set_ptt PTT = true");
          error_check (rig_set_ptt (rig_.data (), RIG_VFO_CURR
                                    , RIG_PTT_RIG_MICDATA == rig_->caps->ptt_type && back_ptt_port_
                                    ? RIG_PTT_ON_DATA : RIG_PTT_ON), tr ("setting PTT on"));
        }
    }
  else
    {
      if (RIG_PTT_NONE != rig_->state.pttport.type.ptt)
        {
          TRACE_CAT ("HamlibTransceiver", "rig_set_ptt PTT = false");
          error_check (rig_set_ptt (rig_.data (), RIG_VFO_CURR, RIG_PTT_OFF), tr ("setting PTT off"));
        }
    }

  update_PTT (on);
}

void HamlibTransceiver::set_conf (char const * item, char const * value)
{
  token_t token = rig_token_lookup (rig_.data (), item);
  if (RIG_CONF_END != token)	// only set if valid for rig model
    {
      error_check (rig_set_conf (rig_.data (), token, value), tr ("setting a configuration item"));
    }
}

QByteArray HamlibTransceiver::get_conf (char const * item)
{
  token_t token = rig_token_lookup (rig_.data (), item);
  QByteArray value {128, '\0'};
  if (RIG_CONF_END != token)	// only get if valid for rig model
    {
      error_check (rig_get_conf (rig_.data (), token, value.data ()), tr ("getting a configuration item"));
    }
  return value;
}

auto HamlibTransceiver::map_mode (rmode_t m) const -> MODE
{
  switch (m)
    {
    case RIG_MODE_AM:
    case RIG_MODE_SAM:
    case RIG_MODE_AMS:
    case RIG_MODE_DSB:
      return AM;

    case RIG_MODE_CW:
      return CW;

    case RIG_MODE_CWR:
      return CW_R;

    case RIG_MODE_USB:
    case RIG_MODE_ECSSUSB:
    case RIG_MODE_SAH:
    case RIG_MODE_FAX:
      return USB;

    case RIG_MODE_LSB:
    case RIG_MODE_ECSSLSB:
    case RIG_MODE_SAL:
      return LSB;

    case RIG_MODE_RTTY:
      return FSK;

    case RIG_MODE_RTTYR:
      return FSK_R;

    case RIG_MODE_PKTLSB:
      return DIG_L;

    case RIG_MODE_PKTUSB:
      return DIG_U;

    case RIG_MODE_FM:
    case RIG_MODE_WFM:
      return FM;

    case RIG_MODE_PKTFM:
      return DIG_FM;

    default:
      return UNK;
    }
}

rmode_t HamlibTransceiver::map_mode (MODE mode) const
{
  switch (mode)
    {
    case AM: return RIG_MODE_AM;
    case CW: return RIG_MODE_CW;
    case CW_R: return RIG_MODE_CWR;
    case USB: return RIG_MODE_USB;
    case LSB: return RIG_MODE_LSB;
    case FSK: return RIG_MODE_RTTY;
    case FSK_R: return RIG_MODE_RTTYR;
    case DIG_L: return RIG_MODE_PKTLSB;
    case DIG_U: return RIG_MODE_PKTUSB;
    case FM: return RIG_MODE_FM;
    case DIG_FM: return RIG_MODE_PKTFM;
    default: break;
    }
  return RIG_MODE_USB;	// quieten compiler grumble
}


code:
// Boost.Function library

//  Copyright Douglas Gregor 2002-2003. Use, modification and
//  distribution is subject to the Boost Software License, Version
//  1.0. (See accompanying file LICENSE_1_0.txt or copy at
//  http://www.boost.org/LICENSE_1_0.txt)

// For more information, see http://www.boost.org

#define BOOST_FUNCTION_NUM_ARGS 5
#include <boost/function/detail/maybe_include.hpp>
#undef BOOST_FUNCTION_NUM_ARGS

license: mit


code:
//---------------------------------------------------------------------------//
// Copyright (c) 2013-2015 Kyle Lutz <kyle.r.lutz@gmail.com>
//
// Distributed under the Boost Software License, Version 1.0
// See accompanying file LICENSE_1_0.txt or copy at
// http://www.boost.org/LICENSE_1_0.txt
//
// See http://boostorg.github.com/compute for more information.
//---------------------------------------------------------------------------//

#ifndef BOOST_COMPUTE_DETAIL_PARAMETER_CACHE_HPP
#define BOOST_COMPUTE_DETAIL_PARAMETER_CACHE_HPP

#include <algorithm>
#include <string>

#include <boost/shared_ptr.hpp>
#include <boost/make_shared.hpp>
#include <boost/noncopyable.hpp>

#include <boost/compute/config.hpp>
#include <boost/compute/device.hpp>
#include <boost/compute/detail/global_static.hpp>
#include <boost/compute/version.hpp>

#ifdef BOOST_COMPUTE_USE_OFFLINE_CACHE
#include <boost/algorithm/string/trim.hpp>
#include <boost/compute/detail/path.hpp>
#include <boost/property_tree/ptree.hpp>
#include <boost/property_tree/json_parser.hpp>
#endif // BOOST_COMPUTE_USE_OFFLINE_CACHE

namespace boost {
namespace compute {
namespace detail {

class parameter_cache : boost::noncopyable
{
public:
    parameter_cache(const device &device)
        : m_dirty(false),
          m_device_name(device.name())
    {
    #ifdef BOOST_COMPUTE_USE_OFFLINE_CACHE
        // get offline cache file name (e.g. /home/user/.boost_compute/tune/device.json)
        m_file_name = make_file_name();

        // load parameters from offline cache file (if it exists)
        if(boost::filesystem::exists(m_file_name)){
            read_from_disk();
        }
    #endif // BOOST_COMPUTE_USE_OFFLINE_CACHE
    }

    ~parameter_cache()
    {
    #ifdef BOOST_COMPUTE_USE_OFFLINE_CACHE
        write_to_disk();
    #endif // BOOST_COMPUTE_USE_OFFLINE_CACHE
    }

    void set(const std::string &object, const std::string &parameter, uint_ value)
    {
        m_cache[std::make_pair(object, parameter)] = value;

        // set the dirty flag to true. this will cause the updated parameters
        // to be stored to disk.
        m_dirty = true;
    }

    uint_ get(const std::string &object, const std::string &parameter, uint_ default_value)
    {
        std::map<std::pair<std::string, std::string>, uint_>::iterator
            iter = m_cache.find(std::make_pair(object, parameter));
        if(iter != m_cache.end()){
            return iter->second;
        }
        else {
            return default_value;
        }
    }

    static boost::shared_ptr<parameter_cache> get_global_cache(const device &device)
    {
        // device name -> parameter cache
        typedef std::map<std::string, boost::shared_ptr<parameter_cache> > cache_map;

        BOOST_COMPUTE_DETAIL_GLOBAL_STATIC(cache_map, caches, ((std::less<std::string>())));

        cache_map::iterator iter = caches.find(device.name());
        if(iter == caches.end()){
            boost::shared_ptr<parameter_cache> cache =
                boost::make_shared<parameter_cache>(device);

            caches.insert(iter, std::make_pair(device.name(), cache));

            return cache;
        }
        else {
            return iter->second;
        }
    }

private:
#ifdef BOOST_COMPUTE_USE_OFFLINE_CACHE
    // returns a string containing a cannoical device name
    static std::string cannonical_device_name(std::string name)
    {
        boost::algorithm::trim(name);
        std::replace(name.begin(), name.end(), ' ', '_');
        std::replace(name.begin(), name.end(), '(', '_');
        std::replace(name.begin(), name.end(), ')', '_');
        return name;
    }

    // returns the boost.compute version string
    static std::string version_string()
    {
        char buf[32];
        std::snprintf(buf, sizeof(buf), "%d.%d.%d", BOOST_COMPUTE_VERSION_MAJOR,
                                                    BOOST_COMPUTE_VERSION_MINOR,
                                                    BOOST_COMPUTE_VERSION_PATCH);
        return buf;
    }

    // returns the file path for the cached parameters
    std::string make_file_name() const
    {
        return detail::parameter_cache_path(true) + cannonical_device_name(m_device_name) + ".json";
    }

    // store current parameters to disk
    void write_to_disk()
    {
        BOOST_ASSERT(!m_file_name.empty());

        if(m_dirty){
            // save current parameters to disk
            boost::property_tree::ptree pt;
            pt.put("header.device", m_device_name);
            pt.put("header.version", version_string());
            typedef std::map<std::pair<std::string, std::string>, uint_> map_type;
            for(map_type::const_iterator iter = m_cache.begin(); iter != m_cache.end(); ++iter){
                const std::pair<std::string, std::string> &key = iter->first;
                pt.add(key.first + "." + key.second, iter->second);
            }
            write_json(m_file_name, pt);

            m_dirty = false;
        }
    }

    // load stored parameters from disk
    void read_from_disk()
    {
        BOOST_ASSERT(!m_file_name.empty());

        m_cache.clear();

        boost::property_tree::ptree pt;
        try {
            read_json(m_file_name, pt);
        }
        catch(boost::property_tree::json_parser::json_parser_error&){
            // no saved cache file, ignore
            return;
        }

        std::string stored_device;
        try {
            stored_device = pt.get<std::string>("header.device");
        }
        catch(boost::property_tree::ptree_bad_path&){
            return;
        }

        std::string stored_version;
        try {
            stored_version = pt.get<std::string>("header.version");
        }
        catch(boost::property_tree::ptree_bad_path&){
            return;
        }

        if(stored_device == m_device_name && stored_version == version_string()){
            typedef boost::property_tree::ptree::const_iterator pt_iter;
            for(pt_iter iter = pt.begin(); iter != pt.end(); ++iter){
                if(iter->first == "header"){
                    // skip header
                    continue;
                }

                boost::property_tree::ptree child_pt = pt.get_child(iter->first);
                for(pt_iter child_iter = child_pt.begin(); child_iter != child_pt.end(); ++child_iter){
                    set(iter->first, child_iter->first, boost::lexical_cast<uint_>(child_iter->second.data()));
                }
            }
        }

        m_dirty = false;
    }
#endif // BOOST_COMPUTE_USE_OFFLINE_CACHE

private:
    bool m_dirty;
    std::string m_device_name;
    std::string m_file_name;
    std::map<std::pair<std::string, std::string>, uint_> m_cache;
};

} // end detail namespace
} // end compute namespace
} // end boost namespace

#endif // BOOST_COMPUTE_DETAIL_PARAMETER_CACHE_HPP


code:
subroutine unpackprop(n1,k,muf,ccur,cxp)

  character ccur*4,cxp*2

  muf=mod(n1,62)
  n1=n1/62

  k=mod(n1,11)
  n1=n1/11

  j=mod(n1,53)
  n1=n1/53
  if(j.eq.0) cxp='*'
  if(j.ge.1 .and. j.le.26) cxp=char(64+j)
  if(j.gt.26) cxp=char(64+j-26)//char(64+j-26)

  j=mod(n1,53)
  n1=n1/53
  if(j.eq.0) ccur(2:2)='*'
  if(j.ge.1 .and. j.le.26) ccur(2:2)=char(64+j)
  if(j.gt.26) ccur(2:3)=char(64+j-26)//char(64+j-26)
  j=n1
  if(j.eq.0) ccur(1:1)='*'
  if(j.ge.1 .and. j.le.26) ccur(1:1)=char(64+j)
  if(j.gt.26) ccur=char(64+j-26)//char(64+j-26)//ccur(2:3)

  return
end subroutine unpackprop


code:
//  Copyright John Maddock 2010, 2012.
//  Copyright Paul A. Bristow 2011, 2012.

//  Use, modification and distribution are subject to the
//  Boost Software License, Version 1.0. (See accompanying file
//  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)

#ifndef BOOST_MATH_CALCULATE_CONSTANTS_CONSTANTS_INCLUDED
#define BOOST_MATH_CALCULATE_CONSTANTS_CONSTANTS_INCLUDED

#include <boost/math/special_functions/trunc.hpp>

namespace boost{ namespace math{ namespace constants{ namespace detail{

template <class T>
template<int N>
inline T constant_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING

   return ldexp(acos(T(0)), 1);

   /*
   // Although this code works well, it's usually more accurate to just call acos
   // and access the number types own representation of PI which is usually calculated
   // at slightly higher precision...

   T result;
   T a = 1;
   T b;
   T A(a);
   T B = 0.5f;
   T D = 0.25f;

   T lim;
   lim = boost::math::tools::epsilon<T>();

   unsigned k = 1;

   do
   {
      result = A + B;
      result = ldexp(result, -2);
      b = sqrt(B);
      a += b;
      a = ldexp(a, -1);
      A = a * a;
      B = A - result;
      B = ldexp(B, 1);
      result = A - B;
      bool neg = boost::math::sign(result) < 0;
      if(neg)
         result = -result;
      if(result <= lim)
         break;
      if(neg)
         result = -result;
      result = ldexp(result, k - 1);
      D -= result;
      ++k;
      lim = ldexp(lim, 1);
   }
   while(true);

   result = B / D;
   return result;
   */
}

template <class T>
template<int N>
inline T constant_two_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return 2 * pi<T, policies::policy<policies::digits2<N> > >();
}

template <class T> // 2 / pi
template<int N>
inline T constant_two_div_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return 2 / pi<T, policies::policy<policies::digits2<N> > >();
}

template <class T>  // sqrt(2/pi)
template <int N>
inline T constant_root_two_div_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt((2 / pi<T, policies::policy<policies::digits2<N> > >()));
}

template <class T>
template<int N>
inline T constant_one_div_two_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return 1 / two_pi<T, policies::policy<policies::digits2<N> > >();
}

template <class T>
template<int N>
inline T constant_root_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(pi<T, policies::policy<policies::digits2<N> > >());
}

template <class T>
template<int N>
inline T constant_root_half_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(pi<T, policies::policy<policies::digits2<N> > >() / 2);
}

template <class T>
template<int N>
inline T constant_root_two_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(two_pi<T, policies::policy<policies::digits2<N> > >());
}

template <class T>
template<int N>
inline T constant_log_root_two_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return log(root_two_pi<T, policies::policy<policies::digits2<N> > >());
}

template <class T>
template<int N>
inline T constant_root_ln_four<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(log(static_cast<T>(4)));
}

template <class T>
template<int N>
inline T constant_e<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   //
   // Although we can clearly calculate this from first principles, this hooks into
   // T's own notion of e, which hopefully will more accurate than one calculated to
   // a few epsilon:
   //
   BOOST_MATH_STD_USING
   return exp(static_cast<T>(1));
}

template <class T>
template<int N>
inline T constant_half<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return static_cast<T>(1) / static_cast<T>(2);
}

template <class T>
template<int M>
inline T constant_euler<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<M>))
{
   BOOST_MATH_STD_USING
   //
   // This is the method described in:
   // "Some New Algorithms for High-Precision Computation of Euler's Constant"
   // Richard P Brent and Edwin M McMillan.
   // Mathematics of Computation, Volume 34, Number 149, Jan 1980, pages 305-312.
   // See equation 17 with p = 2.
   //
   T n = 3 + (M ? (std::min)(M, tools::digits<T>()) : tools::digits<T>()) / 4;
   T lim = M ? ldexp(T(1), 1 - (std::min)(M, tools::digits<T>())) : tools::epsilon<T>();
   T lnn = log(n);
   T term = 1;
   T N = -lnn;
   T D = 1;
   T Hk = 0;
   T one = 1;

   for(unsigned k = 1;; ++k)
   {
      term *= n * n;
      term /= k * k;
      Hk += one / k;
      N += term * (Hk - lnn);
      D += term;

      if(term < D * lim)
         break;
   }
   return N / D;
}

template <class T>
template<int N>
inline T constant_euler_sqr<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
  BOOST_MATH_STD_USING
  return euler<T, policies::policy<policies::digits2<N> > >()
     * euler<T, policies::policy<policies::digits2<N> > >();
}

template <class T>
template<int N>
inline T constant_one_div_euler<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
  BOOST_MATH_STD_USING
  return static_cast<T>(1)
     / euler<T, policies::policy<policies::digits2<N> > >();
}

template <class T>
template<int N>
inline T constant_root_two<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(static_cast<T>(2));
}

template <class T>
template<int N>
inline T constant_root_three<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(static_cast<T>(3));
}

template <class T>
template<int N>
inline T constant_half_root_two<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(static_cast<T>(2)) / 2;
}

template <class T>
template<int N>
inline T constant_ln_two<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   //
   // Although there are good ways to calculate this from scratch, this hooks into
   // T's own notion of log(2) which will hopefully be accurate to the full precision
   // of T:
   //
   BOOST_MATH_STD_USING
   return log(static_cast<T>(2));
}

template <class T>
template<int N>
inline T constant_ln_ten<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return log(static_cast<T>(10));
}

template <class T>
template<int N>
inline T constant_ln_ln_two<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return log(log(static_cast<T>(2)));
}

template <class T>
template<int N>
inline T constant_third<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return static_cast<T>(1) / static_cast<T>(3);
}

template <class T>
template<int N>
inline T constant_twothirds<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return static_cast<T>(2) / static_cast<T>(3);
}

template <class T>
template<int N>
inline T constant_two_thirds<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return static_cast<T>(2) / static_cast<T>(3);
}

template <class T>
template<int N>
inline T constant_three_quarters<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return static_cast<T>(3) / static_cast<T>(4);
}

template <class T>
template<int N>
inline T constant_pi_minus_three<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return pi<T, policies::policy<policies::digits2<N> > >() - static_cast<T>(3);
}

template <class T>
template<int N>
inline T constant_four_minus_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return static_cast<T>(4) - pi<T, policies::policy<policies::digits2<N> > >();
}

//template <class T>
//template<int N>
//inline T constant_pow23_four_minus_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
//{
//   BOOST_MATH_STD_USING
//   return pow(four_minus_pi<T, policies::policy<policies::digits2<N> > >(), static_cast<T>(1.5));
//}

template <class T>
template<int N>
inline T constant_exp_minus_half<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return exp(static_cast<T>(-0.5));
}

// Pi
template <class T>
template<int N>
inline T constant_one_div_root_two<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return static_cast<T>(1) / root_two<T, policies::policy<policies::digits2<N> > >();
}

template <class T>
template<int N>
inline T constant_one_div_root_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return static_cast<T>(1) / root_pi<T, policies::policy<policies::digits2<N> > >();
}

template <class T>
template<int N>
inline T constant_one_div_root_two_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   return static_cast<T>(1) / root_two_pi<T, policies::policy<policies::digits2<N> > >();
}

template <class T>
template<int N>
inline T constant_root_one_div_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(static_cast<T>(1) / pi<T, policies::policy<policies::digits2<N> > >());
}

template <class T>
template<int N>
inline T constant_four_thirds_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >() * static_cast<T>(4) / static_cast<T>(3);
}

template <class T>
template<int N>
inline T constant_half_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >()  / static_cast<T>(2);
}

template <class T>
template<int N>
inline T constant_third_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >()  / static_cast<T>(3);
}

template <class T>
template<int N>
inline T constant_sixth_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >()  / static_cast<T>(6);
}

template <class T>
template<int N>
inline T constant_two_thirds_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >() * static_cast<T>(2) / static_cast<T>(3);
}

template <class T>
template<int N>
inline T constant_three_quarters_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >() * static_cast<T>(3) / static_cast<T>(4);
}

template <class T>
template<int N>
inline T constant_pi_pow_e<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pow(pi<T, policies::policy<policies::digits2<N> > >(), e<T, policies::policy<policies::digits2<N> > >()); //
}

template <class T>
template<int N>
inline T constant_pi_sqr<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >()
   *      pi<T, policies::policy<policies::digits2<N> > >() ; //
}

template <class T>
template<int N>
inline T constant_pi_sqr_div_six<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >()
   *      pi<T, policies::policy<policies::digits2<N> > >()
   / static_cast<T>(6); //
}

template <class T>
template<int N>
inline T constant_pi_cubed<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >()
   *      pi<T, policies::policy<policies::digits2<N> > >()
   *      pi<T, policies::policy<policies::digits2<N> > >()
   ; //
}

template <class T>
template<int N>
inline T constant_cbrt_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pow(pi<T, policies::policy<policies::digits2<N> > >(), static_cast<T>(1)/ static_cast<T>(3));
}

template <class T>
template<int N>
inline T constant_one_div_cbrt_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return static_cast<T>(1)
   / pow(pi<T, policies::policy<policies::digits2<N> > >(), static_cast<T>(1)/ static_cast<T>(3));
}

// Euler's e

template <class T>
template<int N>
inline T constant_e_pow_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pow(e<T, policies::policy<policies::digits2<N> > >(), pi<T, policies::policy<policies::digits2<N> > >()); //
}

template <class T>
template<int N>
inline T constant_root_e<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sqrt(e<T, policies::policy<policies::digits2<N> > >());
}

template <class T>
template<int N>
inline T constant_log10_e<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return log10(e<T, policies::policy<policies::digits2<N> > >());
}

template <class T>
template<int N>
inline T constant_one_div_log10_e<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return  static_cast<T>(1) /
     log10(e<T, policies::policy<policies::digits2<N> > >());
}

// Trigonometric

template <class T>
template<int N>
inline T constant_degree<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return pi<T, policies::policy<policies::digits2<N> > >()
   / static_cast<T>(180)
   ; //
}

template <class T>
template<int N>
inline T constant_radian<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return static_cast<T>(180)
   / pi<T, policies::policy<policies::digits2<N> > >()
   ; //
}

template <class T>
template<int N>
inline T constant_sin_one<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sin(static_cast<T>(1)) ; //
}

template <class T>
template<int N>
inline T constant_cos_one<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return cos(static_cast<T>(1)) ; //
}

template <class T>
template<int N>
inline T constant_sinh_one<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return sinh(static_cast<T>(1)) ; //
}

template <class T>
template<int N>
inline T constant_cosh_one<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return cosh(static_cast<T>(1)) ; //
}

template <class T>
template<int N>
inline T constant_phi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return (static_cast<T>(1) + sqrt(static_cast<T>(5)) )/static_cast<T>(2) ; //
}

template <class T>
template<int N>
inline T constant_ln_phi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   //return  log(phi<T, policies::policy<policies::digits2<N> > >()); // ???
   return log((static_cast<T>(1) + sqrt(static_cast<T>(5)) )/static_cast<T>(2) );
}
template <class T>
template<int N>
inline T constant_one_div_ln_phi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING
   return static_cast<T>(1) /
     log((static_cast<T>(1) + sqrt(static_cast<T>(5)) )/static_cast<T>(2) );
}

// Zeta

template <class T>
template<int N>
inline T constant_zeta_two<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   BOOST_MATH_STD_USING

     return pi<T, policies::policy<policies::digits2<N> > >()
     *  pi<T, policies::policy<policies::digits2<N> > >()
     /static_cast<T>(6);
}

template <class T>
template<int N>
inline T constant_zeta_three<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   // http://mathworld.wolfram.com/AperysConstant.html
   // http://en.wikipedia.org/wiki/Mathematical_constant

   // http://oeis.org/A002117/constant
   //T zeta3("1.20205690315959428539973816151144999076"
   // "4986292340498881792271555341838205786313"
   // "09018645587360933525814619915");

  //"1.202056903159594285399738161511449990, 76498629234049888179227155534183820578631309018645587360933525814619915"  A002117
  // 1.202056903159594285399738161511449990, 76498629234049888179227155534183820578631309018645587360933525814619915780, +00);
  //"1.2020569031595942 double
  // http://www.spaennare.se/SSPROG/ssnum.pdf  // section 11, Algorithm for Apery's constant zeta(3).
  // Programs to Calculate some Mathematical Constants to Large Precision, Document Version 1.50

  // by Stefan Spannare  September 19, 2007
  // zeta(3) = 1/64 * sum
   BOOST_MATH_STD_USING
   T n_fact=static_cast<T>(1); // build n! for n = 0.
   T sum = static_cast<double>(77); // Start with n = 0 case.
   // for n = 0, (77/1) /64 = 1.203125
   //double lim = std::numeric_limits<double>::epsilon();
   T lim = N ? ldexp(T(1), 1 - (std::min)(N, tools::digits<T>())) : tools::epsilon<T>();
   for(unsigned int n = 1; n < 40; ++n)
   { // three to five decimal digits per term, so 40 should be plenty for 100 decimal digits.
      //cout << "n = " << n << endl;
      n_fact *= n; // n!
      T n_fact_p10 = n_fact * n_fact * n_fact * n_fact * n_fact * n_fact * n_fact * n_fact * n_fact * n_fact; // (n!)^10
      T num = ((205 * n * n) + (250 * n) + 77) * n_fact_p10; // 205n^2 + 250n + 77
      // int nn = (2 * n + 1);
      // T d = factorial(nn); // inline factorial.
      T d = 1;
      for(unsigned int i = 1; i <= (n+n + 1); ++i) // (2n + 1)
      {
        d *= i;
      }
      T den = d * d * d * d * d; // [(2n+1)!]^5
      //cout << "den = " << den << endl;
      T term = num/den;
      if (n % 2 != 0)
      { //term *= -1;
        sum -= term;
      }
      else
      {
        sum += term;
      }
      //cout << "term = " << term << endl;
      //cout << "sum/64  = " << sum/64 << endl;
      if(abs(term) < lim)
      {
         break;
      }
   }
   return sum / 64;
}

template <class T>
template<int N>
inline T constant_catalan<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{ // http://oeis.org/A006752/constant
  //T c("0.915965594177219015054603514932384110774"
 //"149374281672134266498119621763019776254769479356512926115106248574");

  // 9.159655941772190150546035149323841107, 74149374281672134266498119621763019776254769479356512926115106248574422619, -01);

   // This is equation (entry) 31 from
   // http://www-2.cs.cmu.edu/~adamchik/articles/catalan/catalan.htm
   // See also http://www.mpfr.org/algorithms.pdf
   BOOST_MATH_STD_USING
   T k_fact = 1;
   T tk_fact = 1;
   T sum = 1;
   T term;
   T lim = N ? ldexp(T(1), 1 - (std::min)(N, tools::digits<T>())) : tools::epsilon<T>();

   for(unsigned k = 1;; ++k)
   {
      k_fact *= k;
      tk_fact *= (2 * k) * (2 * k - 1);
      term = k_fact * k_fact / (tk_fact * (2 * k + 1) * (2 * k + 1));
      sum += term;
      if(term < lim)
      {
         break;
      }
   }
   return boost::math::constants::pi<T, boost::math::policies::policy<> >()
      * log(2 + boost::math::constants::root_three<T, boost::math::policies::policy<> >())
       / 8
      + 3 * sum / 8;
}

namespace khinchin_detail{

template <class T>
T zeta_polynomial_series(T s, T sc, int digits)
{
   BOOST_MATH_STD_USING
   //
   // This is algorithm 3 from:
   //
   // "An Efficient Algorithm for the Riemann Zeta Function", P. Borwein,
   // Canadian Mathematical Society, Conference Proceedings, 2000.
   // See: http://www.cecm.sfu.ca/personal/pborwein/PAPERS/P155.pdf
   //
   BOOST_MATH_STD_USING
   int n = (digits * 19) / 53;
   T sum = 0;
   T two_n = ldexp(T(1), n);
   int ej_sign = 1;
   for(int j = 0; j < n; ++j)
   {
      sum += ej_sign * -two_n / pow(T(j + 1), s);
      ej_sign = -ej_sign;
   }
   T ej_sum = 1;
   T ej_term = 1;
   for(int j = n; j <= 2 * n - 1; ++j)
   {
      sum += ej_sign * (ej_sum - two_n) / pow(T(j + 1), s);
      ej_sign = -ej_sign;
      ej_term *= 2 * n - j;
      ej_term /= j - n + 1;
      ej_sum += ej_term;
   }
   return -sum / (two_n * (1 - pow(T(2), sc)));
}

template <class T>
T khinchin(int digits)
{
   BOOST_MATH_STD_USING
   T sum = 0;
   T term;
   T lim = ldexp(T(1), 1-digits);
   T factor = 0;
   unsigned last_k = 1;
   T num = 1;
   for(unsigned n = 1;; ++n)
   {
      for(unsigned k = last_k; k <= 2 * n - 1; ++k)
      {
         factor += num / k;
         num = -num;
      }
      last_k = 2 * n;
      term = (zeta_polynomial_series(T(2 * n), T(1 - T(2 * n)), digits) - 1) * factor / n;
      sum += term;
      if(term < lim)
         break;
   }
   return exp(sum / boost::math::constants::ln_two<T, boost::math::policies::policy<> >());
}

}

template <class T>
template<int N>
inline T constant_khinchin<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{
   int n = N ? (std::min)(N, tools::digits<T>()) : tools::digits<T>();
   return khinchin_detail::khinchin<T>(n);
}

template <class T>
template<int N>
inline T constant_extreme_value_skewness<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{ // from e_float constants.cpp
  // Mathematica: N[12 Sqrt[6]  Zeta[3]/Pi^3, 1101]
   BOOST_MATH_STD_USING
   T ev(12 * sqrt(static_cast<T>(6)) * zeta_three<T, policies::policy<policies::digits2<N> > >()
    / pi_cubed<T, policies::policy<policies::digits2<N> > >() );

//T ev(
//"1.1395470994046486574927930193898461120875997958365518247216557100852480077060706857071875468869385150"
//"1894272048688553376986765366075828644841024041679714157616857834895702411080704529137366329462558680"
//"2015498788776135705587959418756809080074611906006528647805347822929577145038743873949415294942796280"
//"0895597703063466053535550338267721294164578901640163603544404938283861127819804918174973533694090594"
//"3094963822672055237678432023017824416203652657301470473548274848068762500300316769691474974950757965"
//"8640779777748741897542093874605477776538884083378029488863880220988107155275203245233994097178778984"
//"3488995668362387892097897322246698071290011857605809901090220903955815127463328974447572119951192970"
//"3684453635456559086126406960279692862247058250100678008419431185138019869693206366891639436908462809"
//"9756051372711251054914491837034685476095423926553367264355374652153595857163724698198860485357368964"
//"3807049634423621246870868566707915720704996296083373077647528285782964567312903914752617978405994377"
//"9064157147206717895272199736902453130842229559980076472936976287378945035706933650987259357729800315");

  return ev;
}

namespace detail{
//
// Calculation of the Glaisher constant depends upon calculating the
// derivative of the zeta function at 2, we can then use the relation:
// zeta'(2) = 1/6 pi^2 [euler + ln(2pi)-12ln(A)]
// To get the constant A.
// See equation 45 at http://mathworld.wolfram.com/RiemannZetaFunction.html.
//
// The derivative of the zeta function is computed by direct differentiation
// of the relation:
// (1-2^(1-s))zeta(s) = SUM(n=0, INF){ (-n)^n / (n+1)^s  }
// Which gives us 2 slowly converging but alternating sums to compute,
// for this we use Algorithm 1 from "Convergent Acceleration of Alternating Series",
// Henri Cohen, Fernando Rodriguez Villegas and Don Zagier, Experimental Mathematics 9:1 (1999).
// See http://www.math.utexas.edu/users/villegas/publications/conv-accel.pdf
//
template <class T>
T zeta_series_derivative_2(unsigned digits)
{
   // Derivative of the series part, evaluated at 2:
   BOOST_MATH_STD_USING
   int n = digits * 301 * 13 / 10000;
   boost::math::itrunc((std::numeric_limits<T>::digits10 + 1) * 1.3);
   T d = pow(3 + sqrt(T(8)), n);
   d = (d + 1 / d) / 2;
   T b = -1;
   T c = -d;
   T s = 0;
   for(int k = 0; k < n; ++k)
   {
      T a = -log(T(k+1)) / ((k+1) * (k+1));
      c = b - c;
      s = s + c * a;
      b = (k + n) * (k - n) * b / ((k + T(0.5f)) * (k + 1));
   }
   return s / d;
}

template <class T>
T zeta_series_2(unsigned digits)
{
   // Series part of zeta at 2:
   BOOST_MATH_STD_USING
   int n = digits * 301 * 13 / 10000;
   T d = pow(3 + sqrt(T(8)), n);
   d = (d + 1 / d) / 2;
   T b = -1;
   T c = -d;
   T s = 0;
   for(int k = 0; k < n; ++k)
   {
      T a = T(1) / ((k + 1) * (k + 1));
      c = b - c;
      s = s + c * a;
      b = (k + n) * (k - n) * b / ((k + T(0.5f)) * (k + 1));
   }
   return s / d;
}

template <class T>
inline T zeta_series_lead_2()
{
   // lead part at 2:
   return 2;
}

template <class T>
inline T zeta_series_derivative_lead_2()
{
   // derivative of lead part at 2:
   return -2 * boost::math::constants::ln_two<T>();
}

template <class T>
inline T zeta_derivative_2(unsigned n)
{
   // zeta derivative at 2:
   return zeta_series_derivative_2<T>(n) * zeta_series_lead_2<T>()
      + zeta_series_derivative_lead_2<T>() * zeta_series_2<T>(n);
}

}  // namespace detail

template <class T>
template<int N>
inline T constant_glaisher<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{

   BOOST_MATH_STD_USING
   typedef policies::policy<policies::digits2<N> > forwarding_policy;
   int n = N ? (std::min)(N, tools::digits<T>()) : tools::digits<T>();
   T v = detail::zeta_derivative_2<T>(n);
   v *= 6;
   v /= boost::math::constants::pi<T, forwarding_policy>() * boost::math::constants::pi<T, forwarding_policy>();
   v -= boost::math::constants::euler<T, forwarding_policy>();
   v -= log(2 * boost::math::constants::pi<T, forwarding_policy>());
   v /= -12;
   return exp(v);

 /*
   // from http://mpmath.googlecode.com/svn/data/glaisher.txt
     // 20,000 digits of the Glaisher-Kinkelin constant A = exp(1/2 - zeta'(-1))
     // Computed using A = exp((6 (-zeta'(2))/pi^2 + log 2 pi + gamma)/12)
  // with Euler-Maclaurin summation for zeta'(2).
  T g(
  "1.282427129100622636875342568869791727767688927325001192063740021740406308858826"
  "46112973649195820237439420646120399000748933157791362775280404159072573861727522"
  "14334327143439787335067915257366856907876561146686449997784962754518174312394652"
  "76128213808180219264516851546143919901083573730703504903888123418813674978133050"
  "93770833682222494115874837348064399978830070125567001286994157705432053927585405"
  "81731588155481762970384743250467775147374600031616023046613296342991558095879293"
  "36343887288701988953460725233184702489001091776941712153569193674967261270398013"
  "52652668868978218897401729375840750167472114895288815996668743164513890306962645"
  "59870469543740253099606800842447417554061490189444139386196089129682173528798629"
  "88434220366989900606980888785849587494085307347117090132667567503310523405221054"
  "14176776156308191919997185237047761312315374135304725819814797451761027540834943"
  "14384965234139453373065832325673954957601692256427736926358821692159870775858274"
  "69575162841550648585890834128227556209547002918593263079373376942077522290940187");

  return g;
  */
}

template <class T>
template<int N>
inline T constant_rayleigh_skewness<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{  // From e_float
  // 1100 digits of the Rayleigh distribution skewness
  // Mathematica: N[2 Sqrt[Pi] (Pi - 3)/((4 - Pi)^(3/2)), 1100]

   BOOST_MATH_STD_USING
   T rs(2 * root_pi<T, policies::policy<policies::digits2<N> > >()
      * pi_minus_three<T, policies::policy<policies::digits2<N> > >()
      / pow(four_minus_pi<T, policies::policy<policies::digits2<N> > >(), static_cast<T>(3./2))
      );
 //   6.31110657818937138191899351544227779844042203134719497658094585692926819617473725459905027032537306794400047264,

  //"0.6311106578189371381918993515442277798440422031347194976580945856929268196174737254599050270325373067"
  //"9440004726436754739597525250317640394102954301685809920213808351450851396781817932734836994829371322"
  //"5797376021347531983451654130317032832308462278373358624120822253764532674177325950686466133508511968"
  //"2389168716630349407238090652663422922072397393006683401992961569208109477307776249225072042971818671"
  //"4058887072693437217879039875871765635655476241624825389439481561152126886932506682176611183750503553"
  //"1218982627032068396407180216351425758181396562859085306247387212297187006230007438534686340210168288"
  //"8956816965453815849613622117088096547521391672977226658826566757207615552041767516828171274858145957"
  //"6137539156656005855905288420585194082284972984285863898582313048515484073396332610565441264220790791"
  //"0194897267890422924599776483890102027823328602965235306539844007677157873140562950510028206251529523"
  //"7428049693650605954398446899724157486062545281504433364675815915402937209673727753199567661561209251"
  //"4695589950526053470201635372590001578503476490223746511106018091907936826431407434894024396366284848");  ;
  return rs;
}

template <class T>
template<int N>
inline T constant_rayleigh_kurtosis_excess<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{ // - (6 Pi^2 - 24 Pi + 16)/((Pi - 4)^2)
    // Might provide and calculate this using pi_minus_four.
   BOOST_MATH_STD_USING
   return - (((static_cast<T>(6) * pi<T, policies::policy<policies::digits2<N> > >()
        * pi<T, policies::policy<policies::digits2<N> > >())
   - (static_cast<T>(24) * pi<T, policies::policy<policies::digits2<N> > >()) + static_cast<T>(16) )
   /
   ((pi<T, policies::policy<policies::digits2<N> > >() - static_cast<T>(4))
   * (pi<T, policies::policy<policies::digits2<N> > >() - static_cast<T>(4)))
   );
}

template <class T>
template<int N>
inline T constant_rayleigh_kurtosis<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
{ // 3 - (6 Pi^2 - 24 Pi + 16)/((Pi - 4)^2)
  // Might provide and calculate this using pi_minus_four.
   BOOST_MATH_STD_USING
   return static_cast<T>(3) - (((static_cast<T>(6) * pi<T, policies::policy<policies::digits2<N> > >()
        * pi<T, policies::policy<policies::digits2<N> > >())
   - (static_cast<T>(24) * pi<T, policies::policy<policies::digits2<N> > >()) + static_cast<T>(16) )
   /
   ((pi<T, policies::policy<policies::digits2<N> > >() - static_cast<T>(4))
   * (pi<T, policies::policy<policies::digits2<N> > >() - static_cast<T>(4)))
   );
}

}}}} // namespaces

#endif // BOOST_MATH_CALCULATE_CONSTANTS_CONSTANTS_INCLUDED