text:
// C++ headers
#include <string>
#include <cstring>

// boost headers
#include <boost/regex.hpp>

// Root headers
#include <TH1F.h>

// CMSSW headers
#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/Run.h"
#include "FWCore/Framework/interface/LuminosityBlock.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
#include "FWCore/ParameterSet/interface/Registry.h"
#include "DQMServices/Core/interface/DQMStore.h"
#include "DQMServices/Core/interface/DQMEDHarvester.h"
#include "DQMServices/Core/interface/MonitorElement.h"

struct MEPSet {
  std::string name;
  std::string folder;
};

class PSMonitorClient : public DQMEDHarvester {
public:
  explicit PSMonitorClient(edm::ParameterSet const &);
  ~PSMonitorClient() override = default;

  static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
  static void fillMePSetDescription(edm::ParameterSetDescription &pset);

private:
  static MEPSet getHistoPSet(edm::ParameterSet pset);

  std::string m_dqm_path;

  void dqmEndLuminosityBlock(DQMStore::IBooker &booker,
                             DQMStore::IGetter &getter,
                             edm::LuminosityBlock const &,
                             edm::EventSetup const &) override;
  void dqmEndJob(DQMStore::IBooker &booker, DQMStore::IGetter &getter) override;

  void check(DQMStore::IBooker &booker, DQMStore::IGetter &getter);

  MEPSet psColumnVSlumiPSet;
};

PSMonitorClient::PSMonitorClient(edm::ParameterSet const &config)
    : m_dqm_path(config.getUntrackedParameter<std::string>("dqmPath")),
      psColumnVSlumiPSet(getHistoPSet(config.getParameter<edm::ParameterSet>("me"))) {}

MEPSet PSMonitorClient::getHistoPSet(edm::ParameterSet pset) {
  return MEPSet{
      pset.getParameter<std::string>("name"),
      pset.getParameter<std::string>("folder"),
  };
}

void PSMonitorClient::fillMePSetDescription(edm::ParameterSetDescription &pset) {
  pset.add<std::string>("folder", "HLT/PSMonitoring");
  pset.add<std::string>("name", "psColumnVSlumi");
}

void PSMonitorClient::dqmEndJob(DQMStore::IBooker &booker, DQMStore::IGetter &getter) { check(booker, getter); }

void PSMonitorClient::dqmEndLuminosityBlock(DQMStore::IBooker &booker,
                                            DQMStore::IGetter &getter,
                                            edm::LuminosityBlock const &lumi,
                                            edm::EventSetup const &setup) {
  check(booker, getter);
}

#include "FWCore/MessageLogger/interface/MessageLogger.h"
void PSMonitorClient::check(DQMStore::IBooker &booker, DQMStore::IGetter &getter) {
  std::string folder = psColumnVSlumiPSet.folder;
  std::string name = psColumnVSlumiPSet.name;

  getter.setCurrentFolder(folder);
  MonitorElement *psColumnVSlumi = getter.get(psColumnVSlumiPSet.folder + "/" + psColumnVSlumiPSet.name);
  // if no ME available, return
  if (!psColumnVSlumi) {
    edm::LogWarning("PSMonitorClient") << "no " << psColumnVSlumiPSet.name << " ME is available in "
                                       << psColumnVSlumiPSet.folder << std::endl;
    return;
  }

  /*
  TH2F* h = psColumnVSlumi->getTH2F();
  size_t nbinsX = psColumnVSlumi->getNbinsX();
  size_t nbinsY = psColumnVSlumi->getNbinsY();

  for ( size_t ibinY=1; ibinY < nbinsY; ++ibinY )
    std::cout << h->GetXaxis()->GetBinLabel(ibinY) << std::endl;
    for ( size_t ibinX=1; ibinX< nbinsX; ++ibinX )
      if ( psColumnVSlumi->getBinContent(ibinX) )
	std::cout << "ibinX: " << psColumnVSlumi->getBinContent(ibinX) << std::endl;
  */
}

void PSMonitorClient::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
  // The following says we do not know what parameters are allowed so do no validation
  // Please change this to state exactly what you do use, even if it is no parameters
  edm::ParameterSetDescription desc;
  desc.addUntracked<std::string>("dqmPath", "HLT/PSMonitoring");

  edm::ParameterSetDescription mePSet;
  fillMePSetDescription(mePSet);
  desc.add<edm::ParameterSetDescription>("me", mePSet);

  descriptions.add("psMonitorClient", desc);
}

//define this as a plug-in
#include "FWCore/Framework/interface/MakerMacros.h"
DEFINE_FWK_MODULE(PSMonitorClient);

max_stars_repo_name: eric-moreno/cmssw

max_issues_repo_name: eric-moreno/cmssw

max_forks_repo_name: eric-moreno/cmssw


text:
#include "stdafx.h"
#include "PluginManager.h"
#include <iostream>
#include <bitset>
#include "DeviceContainer.h"
#include "HardwareEventDispatcher.h"
#include <boost/filesystem.hpp>
#include "logger.h"
#include <cassert>
PluginManager::PluginManager(boost::asio::io_service& io)
	: m_plugins()
	, m_deviceContainer(nullptr)
	, m_io(io)
	, m_pluginEventLoop(io, boost::posix_time::millisec(16))
	, m_pluginInfo()
	, m_fatalErrorHandler([]() {})
	, m_virtualDevices()
{
	m_pluginEventLoop.SetEvent([this]() {

		if (!this->TickOnce(16)) {//todo: should store 16 as variable
			m_fatalErrorHandler();
		}
	});

	m_pluginEventLoop.Start();

}

std::vector<boost::filesystem::path> findManifests(const boost::filesystem::path& root) {
	namespace fs = boost::filesystem;

	std::vector<fs::path> paths;

	std::unordered_map<std::string, std::vector<fs::directory_entry>> directoryManifestMap;
	try {
		fs::directory_iterator it(root);
		fs::directory_iterator endit;

		while (it != endit) {
			if (fs::is_directory(*it)) {
				std::vector<fs::directory_entry> dirs;

				std::copy(fs::directory_iterator(*it), fs::directory_iterator(), std::back_inserter(dirs));
				directoryManifestMap[it->path().string()] = dirs;

			}
			++it;
		}
	}
	catch (const boost::filesystem::filesystem_error& ec) {
		BOOST_LOG_SEV(clogger::get(), nsvr_severity_error) << "[PluginDiscovery] " << ec.what();
	}

	for (const auto& dir : directoryManifestMap) {
		for (const auto& entry : dir.second) {
			if (fs::is_regular_file(entry) && Parsing::IsProbablyManifest(entry.path().string())) {
				paths.push_back(entry.path());
			}
		}
	}
	return paths;
}

void PluginManager::Discover()
{
	namespace fs = boost::filesystem;

	std::vector<fs::path> paths = findManifests(fs::current_path() / "plugins");

	BOOST_LOG_SEV(clogger::get(), nsvr_severity_info) << "[PluginManager] " << " Found " << paths.size() << " plugin manifests";
	for (const auto& manifest : paths) {
		BOOST_LOG_SEV(clogger::get(), nsvr_severity_info) <<  "--> " << manifest;
	}

	for (const auto& manifest : paths) {
		BOOST_LOG_SEV(clogger::get(), nsvr_severity_info) << "[PluginManager] " << " Attempting to parse " << manifest.stem().string();
		if (auto config = Parsing::ParseConfig(manifest.string())) {
			std::string stem = manifest.stem().string();
			std::string dllName = stem.substr(0, stem.find_last_of('_'));
			m_pluginInfo[dllName].Descriptor = *config;
			m_pluginInfo[dllName].DllPath = manifest.parent_path().string();
			BOOST_LOG_SEV(clogger::get(), nsvr_severity_info) << "--> Success. Assuming dll name is " << dllName;
		}
		else {
			BOOST_LOG_SEV(clogger::get(), nsvr_severity_info) << "--> Failure";

		}
	}
}

void PluginManager::LoadAll()
{
	for (const auto& nameManifestPair : m_pluginInfo) {
		LoadPlugin(nameManifestPair.second.DllPath, nameManifestPair.first);
	}
}

void PluginManager::UnloadAll()
{
	destroyAll();
}

PluginInstance * PluginManager::MakeVirtualPlugin()
{
	assert(m_deviceContainer);

	auto size = m_plugins.size();
	std::string pluginName = "VirtualDeviceHost" + std::to_string(size);
	auto plugin = std::make_unique<PluginInstance>(m_io, pluginName, size);
	auto dispatcher = std::make_unique<HardwareEventDispatcher>(m_io);
	dispatcher->OnDeviceConnected([name = pluginName, this](nsvr_device_id id, PluginInstance* plugin) {
		m_deviceContainer->AddDevice(id, plugin->apis(),  name, plugin->resources());
	});
	plugin->setDispatcher(std::move(dispatcher));

	m_plugins[pluginName] = std::move(plugin);
	return m_plugins.at(pluginName).get();
}

bool PluginManager::TickOnce(uint64_t dt)
{
	std::vector<std::string> crashedPlugins;
	for (auto& plugin : m_plugins) {
		//if return value is false, we should destroy and reinstantiate the plugin
		//taking care to clean up the DeviceContainer somehow

		if (!plugin.second->tick_once(dt)) {
			crashedPlugins.push_back(plugin.first);
		}

	}

	for (auto& crashedPlugin : crashedPlugins) {
		if (m_plugins[crashedPlugin]->Unload()) {
			m_plugins.erase(crashedPlugin);
			BOOST_LOG_SEV(clogger::get(), nsvr_severity_info) << "[PluginManager] Successfully unloaded misbehaving plugin " << crashedPlugin;
		}
		else {
			BOOST_LOG_SEV(clogger::get(), nsvr_severity_info) << "[PluginManager] Unable to unload misbehaving plugin " << crashedPlugin;
			return false;

		}

	}

	return true;
}

bool PluginManager::Reload(const std::string & name)
{
	//Todo: this must be rewritten if we want this functionality

	/*std::vector<DeviceId<local>> toBeRemoved;
	m_deviceContainer.EachDevice([&name, &toBeRemoved](Device* device) {
		if (device->parentPlugin() == name) {
			toBeRemoved.push_back(device->id());
		}
	});

	for (nsvr_device_id id : toBeRemoved) {
		m_deviceContainer.RemoveDevice(id);
	}

	if (m_plugins[name]->Unload()) {
		return m_plugins[name]->Instantiate();
	}*/

	return false;

}

void PluginManager::OnFatalError(std::function<void()> handler)
{
	m_fatalErrorHandler = handler;
}

std::vector<uint32_t> PluginManager::GetPluginIds() const
{
	std::vector<uint32_t> pluginIds;
	for (const auto& kvp : m_plugins) {
		//Don't bother returning virtual device hosts
		if (kvp.first.find("VirtualDeviceHost") != std::string::npos) {
			continue;
		}
		pluginIds.push_back(kvp.second->GetId());
	}
	return pluginIds;
}

void PluginManager::DrawDiagnostics(uint32_t id, nsvr_diagnostics_ui* ui)
{
	for (auto& kvp : m_plugins) {
		if (kvp.second->GetId() == id) {
			diagnostics_api* api = kvp.second->apis().GetApi<diagnostics_api>();
			if (api) {
				api->submit_updatemenu(ui);
			}

			const auto& vdevices = m_pluginInfo.at(kvp.first).Descriptor.vdevices.devices;
			ui->keyval("Virtual devices available", std::to_string(vdevices.size()).c_str());
			for (const auto& device : vdevices) {
				if (ui->button(device.name.c_str())) {
					m_virtualDevices.push_back(instantiateVirtualDevice(kvp.first, device));
				}
			}

			if (ui->button("Remove all virtual devices")) {
				for (const auto& device : m_virtualDevices) {
					m_deviceContainer->RemoveDevice(device.first, device.second);
				}

				m_virtualDevices.clear();
			}
		}
	}
}

boost::optional<std::pair<std::string, PluginManager::PluginInfo>> PluginManager::GetPluginInfo(uint32_t id)
{
	boost::optional<std::string> foundPlugin;
	for (auto& kvp : m_plugins) {
		if (kvp.second->GetId() == id) {
			foundPlugin = kvp.first;
		}
	}

	if (foundPlugin) {
		return std::make_pair(*foundPlugin, m_pluginInfo[*foundPlugin]);
	}

	return boost::none;
}

void PluginManager::SetDeviceContainer(DeviceContainer* devices)
{
	m_deviceContainer = devices;
}

std::unique_ptr<PluginInstance> PluginManager::linkPlugin(const std::string& filename) {

	//auto dispatcher = std::make_unique<HardwareEventDispatcher>(m_io);
	//dispatcher->OnDeviceConnected([this](nsvr_device_id id, PluginApis& apis, const std::string& pluginName) {
	//	m_deviceContainer.AddDevice(id, apis, m_pluginInfo.at(pluginName).Descriptor.bodygraph, pluginName);
	//});

	//dispatcher->OnDeviceDisconnected([this](nsvr_device_id id) {
	//	m_deviceContainer.RemoveDevice(id);
	//});

	auto instance = std::make_unique<PluginInstance>(m_io, filename, m_plugins.size());

	if (instance->Link()) {
		return instance;
	}
	else {
		return std::unique_ptr<PluginInstance>();
	}

}
bool PluginManager::instantiatePlugin(PluginInstance* plugin)
{
	return plugin->Instantiate();
}

bool PluginManager::configurePlugin(PluginInstance* plugin)
{

	if (!plugin->Configure()) {
		BOOST_LOG_SEV(clogger::get(), nsvr_severity_warning) << "Couldn't configure " << plugin->GetFileName();
		return false;
	}

	return true;
}

PluginInstance* PluginManager::GetHardlightPlugin() {
	return m_plugins["HardlightPlugin"].get();
}

bool PluginManager::LoadPlugin(const std::string& searchDirectory, const std::string& dllName)
{

	auto plugin = linkPlugin((boost::filesystem::path(searchDirectory) / dllName).string());
	if (plugin) {

		m_plugins[dllName] = std::move(plugin);

		auto dispatcher = std::make_unique<HardwareEventDispatcher>(m_io);

		dispatcher->OnDeviceConnected([this, dll = dllName](nsvr_device_id id, PluginInstance* plugin) {
			assert(m_deviceContainer);
			plugin->resources()->bodygraphDescriptor = m_pluginInfo.at(dll).Descriptor.bodygraph;

			m_deviceContainer->AddDevice(id, plugin->apis(), dll, plugin->resources());
		});
		dispatcher->OnDeviceDisconnected([this, dll = dllName](nsvr_device_id id) {
			assert(m_deviceContainer);

			m_deviceContainer->RemoveDevice(DeviceId<local>{id}, dll);
		});

		m_plugins[dllName]->setDispatcher(std::move(dispatcher));

	}
	else {
		return false;
	}

	if (!instantiatePlugin(m_plugins[dllName].get())) {
		return false;
	}

	if (!configurePlugin(m_plugins[dllName].get())) {
		return false;
	}

	return true;
}

void PluginManager::destroyAll()
{
	for (auto& plugin : m_plugins) {
		if (!plugin.second->Unload()) {
			BOOST_LOG_SEV(clogger::get(), nsvr_severity_error) << "[PluginManager] Warning: unable to destroy " << plugin.first;
		}
	}

}

std::pair<DeviceId<local>, std::string> PluginManager::instantiateVirtualDevice(const std::string& plugin, const Parsing::VirtualDeviceDescriptor& device)
{
	const auto& originPlugin = m_pluginInfo.at(plugin);

	auto resources = std::make_unique<PluginInstance::DeviceResources>();
	resources->bodygraphDescriptor = originPlugin.Descriptor.bodygraph;
	//for now we use device ID 0. Probably want to pass in a specific ID.
	resources->deviceDescriptor = DeviceDescriptor{ device.name, 0, device.concept };

	std::vector<DefaultBodygraph::association> assocs;

	std::vector<Node> nodes;
	for (const auto& node : device.nodes) {
		nodes.emplace_back(NodeDescriptor{ node.concept, node.name, node.id });
		for (const auto& region : node.regions) {
			assocs.push_back(DefaultBodygraph::association{ region, node.id });
		}
	}
	resources->discoverer = std::make_unique<DefaultNodeDiscoverer>(nodes);
	resources->id = DeviceId<local>(0);
	resources->waveformHaptics = std::make_unique<DefaultWaveform>();
	resources->bufferedHaptics = std::make_unique<DefaultBuffered>();

	resources->bodygraph = std::make_unique<DefaultBodygraph>(assocs);

	PluginInstance* virtualHost = this->MakeVirtualPlugin();
	virtualHost->addDeviceResources(std::move(resources));

	return std::make_pair(DeviceId<local>{0}, virtualHost->GetFileName());

}

ext: cpp

max_stars_repo_name: cwaldren/HL-engine

max_issues_repo_name: cwaldren/HL-engine

max_forks_repo_name: cwaldren/HL-engine


text:
/**
@file

Copyright John Reid 2008, 2013, 2015

*/

#define BOOST_TEST_MODULE seqan adaptors  // Generate main() in this object file

#ifdef _MSC_VER
# pragma warning( disable : 4244 )
#endif //_MSC_VER

#include "myrrh/seqan_boost_range_adaptors.h"

#include <boost/progress.hpp>
#include <boost/test/unit_test.hpp>
#include <boost/test/unit_test_monitor.hpp>
#include <boost/test/parameterized_test.hpp>
#include <boost/assign/list_of.hpp>
#include <boost/foreach.hpp>
#include <boost/iterator/counting_iterator.hpp>

using namespace boost;
using namespace boost::unit_test;
using namespace boost::assign;
using namespace std;

BOOST_AUTO_TEST_CASE( check_string_set_iterator )
{
    cout << "*********************** check_string_set_iterator()\n";
    typedef seqan::String< seqan::Dna5 >       string_t;            /**< The string type. */
    typedef seqan::StringSet< string_t >       string_set_t;        /**< The string set type. */

    string_set_t string_set;
    seqan::appendValue( string_set, "ACGTGGACTGCGCTAGGCTCGCGC" );
    BOOST_CHECK_EQUAL( 1u, seqan::length( string_set ) );

    BOOST_FOREACH( string_t & s, string_set ) {
        cout << seqan::length( s ) << "\n";
    }

    BOOST_FOREACH( const string_t & s, const_cast< const string_set_t & >( string_set ) ) {
        cout << seqan::length( s ) << "\n";
    }
}

max_stars_repo_name: JohnReid/myrrh

max_issues_repo_name: JohnReid/myrrh

max_forks_repo_name: JohnReid/myrrh


text:
/* -*- Mode:C++; c-file-style:"gnu"; indent-tabs-mode:nil; -*- */
/*
 * Copyright (c) 2013-2019 Regents of the University of California.
 *
 * This file is part of ndn-cxx library (NDN C++ library with eXperimental eXtensions).
 *
 * ndn-cxx library is free software: you can redistribute it and/or modify it under the
 * terms of the GNU Lesser General Public License as published by the Free Software
 * Foundation, either version 3 of the License, or (at your option) any later version.
 *
 * ndn-cxx library is distributed in the hope that it will be useful, but WITHOUT ANY
 * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
 * PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more details.
 *
 * You should have received copies of the GNU General Public License and GNU Lesser
 * General Public License along with ndn-cxx, e.g., in COPYING.md file.  If not, see
 * <http://www.gnu.org/licenses/>.
 *
 * See AUTHORS.md for complete list of ndn-cxx authors and contributors.
 */

#ifndef NDN_DETAIL_OVERLOAD_HPP
#define NDN_DETAIL_OVERLOAD_HPP

#include <boost/version.hpp>

#if BOOST_VERSION >= 106100
#include <boost/hana/functional/overload.hpp>

namespace ndn {
namespace detail {
constexpr boost::hana::make_overload_t overload{};
} // namespace detail
} // namespace ndn

#else
#include <type_traits>

namespace ndn {
namespace detail {
namespace hana_backports {

// The following code is copied from the Boost.Hana library.
// Copyright Louis Dionne 2013-2017
// Distributed under the Boost Software License, Version 1.0.
// (See http://boost.org/LICENSE_1_0.txt)

//! Pick one of several functions to call based on overload resolution.
//!
//! Specifically, `overload(f1, f2, ..., fn)` is a function object such
//! that
//! @code
//!     overload(f1, f2, ..., fn)(x...) == fk(x...)
//! @endcode
//!
//! where `fk` is the function of `f1, ..., fn` that would be called if
//! overload resolution was performed amongst that set of functions only.
//! If more than one function `fk` would be picked by overload resolution,
//! then the call is ambiguous.
#ifndef DOXYGEN
template <typename F, typename ...G>
struct overload_t
    : overload_t<F>::type
    , overload_t<G...>::type
{
  using type = overload_t;
  using overload_t<F>::type::operator();
  using overload_t<G...>::type::operator();

  template <typename F_, typename ...G_>
  constexpr explicit overload_t(F_&& f, G_&& ...g)
    : overload_t<F>::type(static_cast<F_&&>(f))
    , overload_t<G...>::type(static_cast<G_&&>(g)...)
  { }
};

template <typename F>
struct overload_t<F> { using type = F; };

template <typename R, typename ...Args>
struct overload_t<R(*)(Args...)> {
  using type = overload_t;
  R (*fptr_)(Args...);

  explicit constexpr overload_t(R (*fp)(Args...))
    : fptr_(fp)
  { }

  constexpr R operator()(Args ...args) const
  { return fptr_(static_cast<Args&&>(args)...); }
};

struct make_overload_t {
  template <typename ...F,
    typename Overload = typename overload_t<
      typename std::decay<F>::type...
    >::type
  >
  constexpr Overload operator()(F&& ...f) const {
    return Overload(static_cast<F&&>(f)...);
  }
};
#endif // DOXYGEN

} // namespace hana_backports

constexpr hana_backports::make_overload_t overload{};

} // namespace detail
} // namespace ndn

#endif // BOOST_VERSION >= 106100
#endif // NDN_DETAIL_OVERLOAD_HPP

ext: hpp


text:
/**************************************************************
 *
 * Licensed to the Apache Software Foundation (ASF) under one
 * or more contributor license agreements.  See the NOTICE file
 * distributed with this work for additional information
 * regarding copyright ownership.  The ASF licenses this file
 * to you under the Apache License, Version 2.0 (the
 * "License"); you may not use this file except in compliance
 * with the License.  You may obtain a copy of the License at
 *
 *   http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing,
 * software distributed under the License is distributed on an
 * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
 * KIND, either express or implied.  See the License for the
 * specific language governing permissions and limitations
 * under the License.
 *
 *************************************************************/

// MARKER(update_precomp.py): autogen include statement, do not remove
#include "precompiled_basegfx.hxx"
// autogenerated file with codegen.pl

#include "preextstl.h"
#include "gtest/gtest.h"
#include "postextstl.h"

#include <basegfx/tools/b2dclipstate.hxx>
#include <basegfx/range/b2dpolyrange.hxx>
#include <basegfx/polygon/b2dpolygon.hxx>
#include <basegfx/polygon/b2dpolygontools.hxx>
#include <basegfx/polygon/b2dpolypolygontools.hxx>
#include <basegfx/polygon/b2dpolypolygon.hxx>
#include <basegfx/numeric/ftools.hxx>

#include <boost/bind.hpp>

using namespace ::basegfx;

// FIXME:
#define RUN_OLD_FAILING_TESTS 0

namespace basegfx2d
{

class clipstate : public ::testing::Test
{
protected:
    tools::B2DClipState aUnion1;
    tools::B2DClipState aUnion2;
    tools::B2DClipState aIntersect;
    tools::B2DClipState aXor;
    tools::B2DClipState aSubtract;

public:
    virtual void SetUp()
    {
        B2DRange aCenter(100, 100, -100, -100);
        B2DRange aNorth(-10, -110, 10, -90);
        B2DRange aWest(-110, -10, -90, 10);
        B2DRange aSouth(-10, 110, 10, 90);
        B2DRange aEast(110, -10, 90, 10);

        aUnion1.unionRange(aCenter);
        aUnion1.unionRange(aNorth);
        aUnion1.unionRange(aWest);
        aUnion1.unionRange(aSouth);
        aUnion1.unionRange(aEast);

        aUnion2.makeNull();
        aUnion2.unionRange(aCenter);
        aUnion2.unionRange(aNorth);
        aUnion2.unionRange(aWest);
        aUnion2.unionRange(aSouth);
        aUnion2.unionRange(aEast);

        aIntersect.intersectRange(aCenter);
        aIntersect.intersectRange(aNorth);
        aIntersect.intersectRange(aWest);
        aIntersect.intersectRange(aSouth);
        aIntersect.intersectRange(aEast);

        aXor.makeNull();
        aXor.xorRange(aCenter);
        aXor.xorRange(aNorth);
        aXor.xorRange(aWest);
        aXor.xorRange(aSouth);
        aXor.xorRange(aEast);

        aSubtract.intersectRange(aCenter);
        aSubtract.subtractRange(aNorth);
        aSubtract.subtractRange(aWest);
        aSubtract.subtractRange(aSouth);
        aSubtract.subtractRange(aEast);
    }

    virtual void TearDown()
    {}

    void verifyPoly(const char* sName, const char* sSvg, const tools::B2DClipState& toTest)
    {
#if defined(VERBOSE)
        fprintf(stderr, "%s - svg:d=\"%s\"\n",
                sName, rtl::OUStringToOString(
                    basegfx::tools::exportToSvgD(toTest.getClipPoly(), true, true, false),
                    RTL_TEXTENCODING_UTF8).getStr() );
#endif

        B2DPolyPolygon aTmp1;
        ASSERT_TRUE(tools::importFromSvgD(
                               aTmp1, rtl::OUString::createFromAscii(sSvg), false, 0)) << sName;

        const rtl::OUString aSvg=
            tools::exportToSvgD(toTest.getClipPoly(), true, true, false);
        B2DPolyPolygon aTmp2;
        ASSERT_TRUE(tools::importFromSvgD(aTmp2, aSvg, false, 0)) << sName;

        ASSERT_TRUE(aTmp2 == aTmp1) << sName;
    }
};

TEST_F(clipstate, verifySimpleRange)
{
    const char* unionSvg="m100 10v90h-90v10h-20v-10h-90v-90h-10v-20h10v-90h90v-10h20v10h90v90h10v20z";
    const char* intersectSvg="m-100 10v-20h10v20zm80 90v-10h20v10zm-20-190v-10h20v10zm80 100v-20h10v20z";
    const char* xorSvg="m-100 10h10v-20h-10zm90 110h20v-10h-20zm0-180h20v-10h-20zm100 110h10v-20h-10zm10 20v90h-90v10h-20v-10h-90v-90h-10v-20h10v-90h90v-10h20v10h90v90h10v20z";
    const char* subtractSvg="m-90 10v-20h-10v-90h90v10h20v-10h90v90h-10v20h10v90h-90v-10h-20v10h-90v-90z";

    ASSERT_TRUE(aUnion1.isCleared()) << "cleared clip stays empty under union operation";
    verifyPoly("union", unionSvg, aUnion2);
#if RUN_OLD_FAILING_TESTS
    verifyPoly("intersect", intersectSvg, aIntersect);
    verifyPoly("xor", xorSvg, aXor);
#endif
    verifyPoly("subtract", subtractSvg, aSubtract);
}

TEST_F(clipstate, verifyMixedClips)
{
    tools::B2DClipState aMixedClip;

    const char* unionSvg="m100 10v90h-90v10h-20v-10h-90v-90h-10v-20h10v-90h90v-10h20v10h90v90h10v20z";

    B2DPolyPolygon aTmp1;
    tools::importFromSvgD(
        aTmp1, rtl::OUString::createFromAscii(unionSvg), false, 0);

    aMixedClip.intersectPolyPolygon(aTmp1);
    aMixedClip.subtractRange(B2DRange(-20,-150,20,0));
    aMixedClip.subtractRange(B2DRange(-150,-20,0,20));
    aMixedClip.xorRange(B2DRange(-150,-150,150,150));

    const char* mixedClipSvg="m0 0v20h-100v80h90v10h20v-10h90v-90h10v-20h-10v-90h-80v100zm-40-20v-80h-80v80zm-50 170v-300h300v300z";
#if RUN_OLD_FAILING_TESTS
    verifyPoly("mixed clip", mixedClipSvg, aMixedClip);
#endif
}

// -----------------------------------------------------------------------------

} // namespace basegfx2d

ext: cxx

max_stars_repo_name: Grosskopf/openoffice

max_issues_repo_name: Grosskopf/openoffice

max_forks_repo_name: Grosskopf/openoffice


text: #include <boost/test/utils/iterator/input_iterator_facade.hpp>

max_stars_repo_name: miathedev/BoostForArduino

max_issues_repo_name: miathedev/BoostForArduino

max_forks_repo_name: miathedev/BoostForArduino


text:
#include <iostream>
#include <vector>
#include <math.h>
#include <string>

#include "Input_Reader.h"
#include "Materials.h"
#include "Fem_Quadrature.h"
#include "Quadrule_New.h"
#include "Cell_Data.h"
#include "Angular_Quadrature.h"
#include "Time_Data.h"
#include "Temperature_Data.h"

#include <Eigen/Dense>
#include "Diffusion_Matrix_Creator_Grey.h"
#include "Sweep_Matrix_Creator_Grey.h"

#include "Dark_Arts_Exception.h"

/**
  Goal of this unit test is to check source moment formation, and reaction matrix for SLXS Lobatto scheme
  -This is a MMS problem with a spatially varying sigma_a, constant cv , zero sig_s
*/

int main(int argc, char** argv)
{
  int val = 0;
  const double tol = 1.0E-6;

  Input_Reader input_reader;
  try
  {
    input_reader.read_xml(argv[1]);
  }
  catch(const Dark_Arts_Exception& da_exception )
  {
    da_exception.message() ;
    val = -1;
  }

  /// Initialize a Quadrule object to be able to get all of the quadrature we need
  Quadrule_New quad_fun;
  Fem_Quadrature fem_quadrature( input_reader , quad_fun);
  Cell_Data cell_data( input_reader );
  Angular_Quadrature angular_quadrature( input_reader , quad_fun );

  const int n_cell = cell_data.get_total_number_of_cells();
  const int n_dfem_p = fem_quadrature.get_number_of_interpolation_points();

  /// Create a Materials object that contains all opacity, heat capacity, and source objects

  const int n_l_mom = angular_quadrature.get_number_of_leg_moments();
  Materials materials( input_reader, fem_quadrature , cell_data, angular_quadrature );
    Time_Data time_data(input_reader);
    Temperature_Data t_old(fem_quadrature, input_reader, cell_data);
    Temperature_Data t_star(fem_quadrature, input_reader, cell_data);
    Intensity_Data i_old(cell_data, angular_quadrature, fem_quadrature, materials,input_reader);
    K_Temperature kt(n_cell , time_data.get_number_of_stages(), fem_quadrature);
    K_Intensity ki(n_cell ,  time_data.get_number_of_stages(),fem_quadrature, angular_quadrature);
    t_star = t_old;
      const int n_stages = time_data.get_number_of_stages();
    std::shared_ptr<V_Sweep_Matrix_Creator> sweep_matrix;
    sweep_matrix = std::make_shared<Sweep_Matrix_Creator_Grey>(fem_quadrature, materials, n_stages, angular_quadrature.get_sum_w() ,
      n_l_mom, t_old, i_old, kt, ki,t_star) ;

    Diffusion_Matrix_Creator_Grey diffusion_matrix(fem_quadrature, materials, angular_quadrature, t_old, n_cell, input_reader);

    Eigen::MatrixXd sweep_sig_t = Eigen::MatrixXd(n_dfem_p,n_dfem_p);
    Eigen::MatrixXd sweep_sig_s = Eigen::MatrixXd(n_dfem_p,n_dfem_p);
    Eigen::MatrixXd sweep_sig_a;

    Eigen::MatrixXd dsa_sig_s = Eigen::MatrixXd::Zero(n_dfem_p,n_dfem_p);
    Eigen::MatrixXd dsa_sig_a = Eigen::MatrixXd::Zero(n_dfem_p,n_dfem_p);

    const double dt = time_data.get_dt_max();
    const double time_stage = dt + time_data.get_t_start();
    std::vector<double> rk_a_of_stage_i(n_stages,0.);
    rk_a_of_stage_i[0] = 1.;

  try{
    sweep_matrix->set_time_data(dt,time_stage,rk_a_of_stage_i,0);
    diffusion_matrix.set_time_data( dt, time_stage, 1. );

    for(int cell = 0 ; cell < n_cell ; cell++)
    {
      sweep_matrix->update_cell_dependencies(cell);
      sweep_matrix->update_group_dependencies(0);

      std::cout << "Dies here"<<std::endl;
      sweep_matrix->get_r_sig_t(sweep_sig_t);
      sweep_matrix->get_r_sig_s(sweep_sig_s,0);

      sweep_sig_a = sweep_sig_t - sweep_sig_s;

      std::cout << "Here"<<std::endl;

      diffusion_matrix.set_cell_group_information( cell, 0, cell_data.get_cell_width(cell) );
      diffusion_matrix.calculate_pseudo_r_sig_a_and_pseudo_r_sig_s(dsa_sig_a,dsa_sig_s);

      std::cout << "Cell: " << cell << std::endl;
      std::cout << "Transport sweep matrices. \n Sig_s:\n" << sweep_sig_s << "\n Sig_a: \n" << sweep_sig_a << std::endl;

      std::cout << "Diffusion matrices. \n Sig_s:\n" << dsa_sig_s << "\n Sig_a: \n" << dsa_sig_a << std::endl;

      for(int i = 0;i < n_dfem_p ; i++)
      {
        for(int j=0 ; j < n_dfem_p ; j++)
        {
          if( fabs(sweep_sig_a(i,j) - dsa_sig_a(i,j) ) > tol)
            throw Dark_Arts_Exception(MIP , "Difference in R_sig_a");

          if( fabs(sweep_sig_s(i,j) - dsa_sig_s(i,j) ) > tol)
            throw Dark_Arts_Exception(MIP , "Difference in R_sig_s");
        }
      }

    }

  }
  catch(const Dark_Arts_Exception& da)
  {
    val = -1;
    da.testing_message();
  }

  return val;
}

ext: cc


text:
/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */

/*
 Copyright (C) 2014 Peter Caspers

 This file is part of QuantLib, a free-software/open-source library
 for financial quantitative analysts and developers - http://quantlib.org/

 QuantLib is free software: you can redistribute it and/or modify it
 under the terms of the QuantLib license.  You should have received a
 copy of the license along with this program; if not, please email
 <quantlib-dev@lists.sf.net>. The license is also available online at
 <http://quantlib.org/license.shtml>.

 This program is distributed in the hope that it will be useful, but WITHOUT
 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
 FOR A PARTICULAR PURPOSE.  See the license for more details.
*/

#include <ql/experimental/volatility/noarbsabr.hpp>

#include <ql/math/solvers1d/brent.hpp>
#include <ql/math/modifiedbessel.hpp>

#include <boost/math/special_functions/gamma.hpp>
#include <boost/lambda/lambda.hpp>
#include <boost/lambda/bind.hpp>
#include <boost/assign/std/vector.hpp>
#include <boost/functional/hash.hpp>

namespace QuantLib {

class NoArbSabrModel::integrand {
    const NoArbSabrModel* model;
    Real strike;
  public:
    integrand(const NoArbSabrModel* model, Real strike)
    : model(model), strike(strike) {}
    Real operator()(Real f) const {
        return std::max(f - strike, 0.0) * model->p(f);
    }
};

class NoArbSabrModel::p_integrand {
    const NoArbSabrModel* model;
  public:
    explicit p_integrand(const NoArbSabrModel* model)
    : model(model) {}
    Real operator()(Real f) const {
        return model->p(f);
    }
};

NoArbSabrModel::NoArbSabrModel(const Real expiryTime, const Real forward,
                               const Real alpha, const Real beta, const Real nu,
                               const Real rho)
    : expiryTime_(expiryTime), externalForward_(forward), alpha_(alpha),
      beta_(beta), nu_(nu), rho_(rho), forward_(forward),
      numericalForward_(forward) {

    QL_REQUIRE(expiryTime > 0.0 && expiryTime <= detail::NoArbSabrModel::expiryTime_max,
               "expiryTime (" << expiryTime << ") out of bounds");
    QL_REQUIRE(forward > 0.0, "forward (" << forward << ") must be positive");
    QL_REQUIRE(beta >= detail::NoArbSabrModel::beta_min && beta <= detail::NoArbSabrModel::beta_max,
               "beta (" << beta << ") out of bounds");
    Real sigmaI = alpha * std::pow(forward, beta - 1.0);
    QL_REQUIRE(sigmaI >= detail::NoArbSabrModel::sigmaI_min &&
                   sigmaI <= detail::NoArbSabrModel::sigmaI_max,
               "sigmaI = alpha*forward^(beta-1.0) ("
                   << sigmaI << ") out of bounds, alpha=" << alpha
                   << " beta=" << beta << " forward=" << forward);
    QL_REQUIRE(nu >= detail::NoArbSabrModel::nu_min && nu <= detail::NoArbSabrModel::nu_max,
               "nu (" << nu << ") out of bounds");
    QL_REQUIRE(rho >= detail::NoArbSabrModel::rho_min && rho <= detail::NoArbSabrModel::rho_max,
               "rho (" << rho << ") out of bounds");

    // determine a region sufficient for integration in the normal case

    fmin_ = fmax_ = forward_;
    for (Real tmp = p(fmax_);
         tmp > std::max(detail::NoArbSabrModel::i_accuracy / std::max(1.0, fmax_ - fmin_),
                        detail::NoArbSabrModel::density_threshold);
         tmp = p(fmax_)) {
        fmax_ *= 2.0;
    }
    for (Real tmp = p(fmin_);
         tmp > std::max(detail::NoArbSabrModel::i_accuracy / std::max(1.0, fmax_ - fmin_),
                        detail::NoArbSabrModel::density_threshold);
         tmp = p(fmin_)) {
        fmin_ *= 0.5;
    }
    fmin_ = std::max(detail::NoArbSabrModel::strike_min, fmin_);

    QL_REQUIRE(fmax_ > fmin_, "could not find a reasonable integration domain");

    integrator_ =
        ext::make_shared<GaussLobattoIntegral>(
            detail::NoArbSabrModel::i_max_iterations, detail::NoArbSabrModel::i_accuracy);

    detail::D0Interpolator d0(forward_, expiryTime_, alpha_, beta_, nu_, rho_);
    absProb_ = d0();

    try {
        Brent b;
        Real start = std::sqrt(externalForward_ - detail::NoArbSabrModel::strike_min);
        Real tmp =
            b.solve(boost::lambda::bind(&NoArbSabrModel::forwardError, this,
                                        boost::lambda::_1),
                    detail::NoArbSabrModel::forward_accuracy, start,
                    std::min(detail::NoArbSabrModel::forward_search_step, start / 2.0));
        forward_ = tmp * tmp + detail::NoArbSabrModel::strike_min;
    } catch (Error&) {
        // fall back to unadjusted forward
        forward_ = externalForward_;
    }

    Real d = forwardError(std::sqrt(forward_ - detail::NoArbSabrModel::strike_min));
    numericalForward_ = d + externalForward_;
}

Real NoArbSabrModel::optionPrice(const Real strike) const {
    if (p(std::max(forward_, strike)) < detail::NoArbSabrModel::density_threshold)
        return 0.0;
    return (1.0 - absProb_) *
        ((*integrator_)(integrand(this, strike),
                        strike, std::max(fmax_, 2.0 * strike)) /
            numericalIntegralOverP_);
}

Real NoArbSabrModel::digitalOptionPrice(const Real strike) const {
    if (strike < QL_MIN_POSITIVE_REAL)
        return 1.0;
    if (p(std::max(forward_, strike)) < detail::NoArbSabrModel::density_threshold)
        return 0.0;
    return (1.0 - absProb_)
        * ((*integrator_)(p_integrand(this),
                          strike, std::max(fmax_, 2.0 * strike)) /
           numericalIntegralOverP_);
}

Real NoArbSabrModel::forwardError(const Real forward) const {
    forward_ = forward * forward + detail::NoArbSabrModel::strike_min;
    numericalIntegralOverP_ = (*integrator_)(p_integrand(this),
                                             fmin_, fmax_);
    return optionPrice(0.0) - externalForward_;
}

Real NoArbSabrModel::p(const Real f) const {

    if (f < detail::NoArbSabrModel::density_lower_bound ||
        forward_ < detail::NoArbSabrModel::density_lower_bound)
        return 0.0;

    Real fOmB = std::pow(f, 1.0 - beta_);
    Real FOmB = std::pow(forward_, 1.0 - beta_);

    Real zf = fOmB / (alpha_ * (1.0 - beta_));
    Real zF = FOmB / (alpha_ * (1.0 - beta_));
    Real z = zF - zf;

    // Real JzF = std::sqrt(1.0 - 2.0 * rho_ * nu_ * zF + nu_ * nu_ * zF * zF);
    Real Jmzf = std::sqrt(1.0 + 2.0 * rho_ * nu_ * zf + nu_ * nu_ * zf * zf);
    Real Jz = std::sqrt(1.0 - 2.0 * rho_ * nu_ * z + nu_ * nu_ * z * z);

    Real xz = std::log((Jz - rho_ + nu_ * z) / (1.0 - rho_)) / nu_;
    Real Bp_B = beta_ / FOmB;
    // Real Bpp_B = beta_ * (2.0 * beta_ - 1.0) / (FOmB * FOmB);
    Real kappa1 = 0.125 * nu_ * nu_ * (2.0 - 3.0 * rho_ * rho_) -
                  0.25 * rho_ * nu_ * alpha_ * Bp_B;
    // Real kappa2 = alpha_ * alpha_ * (0.25 * Bpp_B - 0.375 * Bp_B * Bp_B);
    Real gamma = 1.0 / (2.0 * (1.0 - beta_));
    Real sqrtOmR = std::sqrt(1.0 - rho_ * rho_);
    Real h = 0.5 * beta_ * rho_ / ((1.0 - beta_) * Jmzf * Jmzf) *
             (nu_ * zf * std::log(zf * Jz / zF) +
              (1 + rho_ * nu_ * zf) / sqrtOmR *
                  (std::atan((nu_ * z - rho_) / sqrtOmR) +
                   std::atan(rho_ / sqrtOmR)));

    Real res =
        std::pow(Jz, -1.5) / (alpha_ * std::pow(f, beta_) * expiryTime_) *
        std::pow(zf, 1.0 - gamma) * std::pow(zF, gamma) *
        std::exp(-(xz * xz) / (2.0 * expiryTime_) +
                 (h + kappa1 * expiryTime_)) *
        modifiedBesselFunction_i_exponentiallyWeighted(gamma,
                                                       zF * zf / expiryTime_);
    return res;
}

namespace detail {

using namespace boost::assign;

D0Interpolator::D0Interpolator(const Real forward, const Real expiryTime,
                               const Real alpha, const Real beta, const Real nu,
                               const Real rho)
    : forward_(forward), expiryTime_(expiryTime), alpha_(alpha), beta_(beta),
      nu_(nu), rho_(rho), gamma_(1.0 / (2.0 * (1.0 - beta_))) {

    sigmaI_ = alpha_ * std::pow(forward_, beta_ - 1.0);

    tauG_ += 0.25, 0.5, 0.75, 1.0, 1.25, 1.5, 1.75, 2.0, 2.25, 2.5, 2.75, 3.0,
        3.25, 3.5, 3.75, 4.0, 4.25, 4.5, 4.75, 5.0, 5.25, 5.5, 5.75, 6.0, 6.25,
        6.5, 6.75, 7.0, 7.25, 7.5, 7.75, 8.0, 8.25, 8.5, 8.75, 9.0, 9.25, 9.5,
        9.75, 10.0, 10.25, 10.5, 10.75, 11.0, 11.25, 11.5, 11.75, 12.0, 12.25,
        12.5, 12.75, 13.0, 13.25, 13.5, 13.75, 14.0, 14.25, 14.5, 14.75, 15.0,
        15.25, 15.5, 15.75, 16.0, 16.25, 16.5, 16.75, 17.0, 17.25, 17.5, 17.75,
        18.0, 18.25, 18.5, 18.75, 19.0, 19.25, 19.5, 19.75, 20.0, 20.25, 20.5,
        20.75, 21.0, 21.25, 21.5, 21.75, 22.0, 22.25, 22.5, 22.75, 23.0, 23.25,
        23.5, 23.75, 24.0, 24.25, 24.5, 24.75, 25.0, 25.25, 25.5, 25.75, 26.0,
        26.25, 26.5, 26.75, 27.0, 27.25, 27.5, 27.75, 28.0, 28.25, 28.5, 28.75,
        29.0, 29.25, 29.5, 29.75, 30.0;

    sigmaIG_ += 1.0, 0.8, 0.7, 0.6, 0.5, 0.45, 0.4, 0.35, 0.3, 0.27, 0.24, 0.21,
        0.18, 0.15, 0.125, 0.1, 0.075, 0.05;

    rhoG_ += 0.75, 0.50, 0.25, 0.00, -0.25, -0.50, -0.75;

    nuG_ += 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8;

    betaG_ += 0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9;
}

Real D0Interpolator::operator()() const {

    // we do not need to check the indices here, because this is already
    // done in the NoArbSabr constructor

    Size tauInd = std::upper_bound(tauG_.begin(), tauG_.end(), expiryTime_) -
                                   tauG_.begin();
    if (tauInd == tauG_.size())
        --tauInd; // tau at upper bound
    Real expiryTimeTmp = expiryTime_;
    if (tauInd == 0) {
        ++tauInd;
        expiryTimeTmp = tauG_.front();
    }
    Real tauL = (expiryTimeTmp - tauG_[tauInd - 1]) /
                (tauG_[tauInd] - tauG_[tauInd - 1]);

    int sigmaIInd =
        sigmaIG_.size() -
        (std::upper_bound(sigmaIG_.rbegin(), sigmaIG_.rend(), sigmaI_) -
         sigmaIG_.rbegin());
    if (sigmaIInd == 0)
        ++sigmaIInd; // sigmaI at upper bound
    Real sigmaIL = (sigmaI_ - sigmaIG_[sigmaIInd - 1]) /
                   (sigmaIG_[sigmaIInd] - sigmaIG_[sigmaIInd - 1]);

    int rhoInd =
        rhoG_.size() -
        (std::upper_bound(rhoG_.rbegin(), rhoG_.rend(), rho_) - rhoG_.rbegin());
    if (rhoInd == 0) {
        rhoInd++;
    }
    if (rhoInd == static_cast<int>(rhoG_.size())) {
        rhoInd--;
    }
    Real rhoL =
        (rho_ - rhoG_[rhoInd - 1]) / (rhoG_[rhoInd] - rhoG_[rhoInd - 1]);

    // for nu = 0 we know phi = 0.5*z_F^2
    Size nuInd = std::upper_bound(nuG_.begin(), nuG_.end(), nu_) - nuG_.begin();
    if (nuInd == nuG_.size())
        --nuInd; // nu at upper bound
    Real tmpNuG = nuInd > 0 ? nuG_[nuInd - 1] : 0.0;
    Real nuL = (nu_ - tmpNuG) / (nuG_[nuInd] - tmpNuG);

    // for beta = 1 we know phi = 0.0
    Size betaInd =
        std::upper_bound(betaG_.begin(), betaG_.end(), beta_) - betaG_.begin();
    Real tmpBetaG;
    if (betaInd == betaG_.size())
        tmpBetaG = 1.0;
    else
        tmpBetaG = betaG_[betaInd];
    Real betaL =
        (beta_ - betaG_[betaInd - 1]) / (tmpBetaG - betaG_[betaInd - 1]);

    Real phiRes = 0.0;
    for (int iTau = -1; iTau <= 0; ++iTau) {
        for (int iSigma = -1; iSigma <= 0; ++iSigma) {
            for (int iRho = -1; iRho <= 0; ++iRho) {
                for (int iNu = -1; iNu <= 0; ++iNu) {
                    for (int iBeta = -1; iBeta <= 0; ++iBeta) {
                        Real phiTmp;
                        if (iNu == -1 && nuInd == 0) {
                            phiTmp =
                                0.5 /
                                (sigmaI_ * sigmaI_ * (1.0 - beta_) *
                                 (1.0 - beta_)); // this is 0.5*z_F^2, see above
                        } else {
                            if (iBeta == 0 && betaInd == betaG_.size()) {
                                phiTmp =
                                    phi(detail::NoArbSabrModel::tiny_prob);
                            } else {
                                int ind = (tauInd + iTau +
                                           (sigmaIInd + iSigma +
                                            (rhoInd + iRho +
                                             (nuInd + iNu + ((betaInd + iBeta) *
                                                             nuG_.size())) *
                                                 rhoG_.size()) *
                                                sigmaIG_.size()) *
                                               tauG_.size());
                                QL_REQUIRE(ind >= 0 && ind < 1209600,
                                           "absorption matrix index ("
                                               << ind << ") invalid");
                                phiTmp = phi((Real)sabrabsprob[ind] /
                                             detail::NoArbSabrModel::nsim);
                            }
                        }
                        phiRes += phiTmp * (iTau == -1 ? (1.0 - tauL) : tauL) *
                                  (iSigma == -1 ? (1.0 - sigmaIL) : sigmaIL) *
                                  (iRho == -1 ? (1.0 - rhoL) : rhoL) *
                                  (iNu == -1 ? (1.0 - nuL) : nuL) *
                                  (iBeta == -1 ? (1.0 - betaL) : betaL);
                    }
                }
            }
        }
    }
    return d0(phiRes);
}

Real D0Interpolator::phi(const Real d0) const {
    if (d0 < 1e-14)
        return detail::NoArbSabrModel::phiByTau_cutoff * expiryTime_;
    return boost::math::gamma_q_inv(gamma_, d0) * expiryTime_;
}

Real D0Interpolator::d0(const Real phi) const {
    return boost::math::gamma_q(gamma_, std::max(0.0, phi / expiryTime_));
}

} // namespace detail

} // namespace QuantLib

max_stars_repo_name: tlapfai/My-Quantlib

max_issues_repo_name: tlapfai/My-Quantlib

max_forks_repo_name: tlapfai/My-Quantlib


text:
#include "PluginProcessor.h"
#include "PluginEditor.h"
#include "AboutDialog.h"
#include <algorithm>
#include <dlib/image_processing/frontal_face_detector.h>
#include <dlib/image_io.h>
#include <dlib/image_processing.h>

#define WAHWTH_USE_IMAGE_PROCESSING_THREAD 1

bool kShowFace = true;

constexpr int kMinWidth = 400;
constexpr int kMinHeight = 300;
constexpr int kMaxWidth = 1200;
constexpr int kMaxHeight = 900;
constexpr int kDefaultWidth = 600;
constexpr int kDefaultHeight = 450;
constexpr int kProcessImageWidth = kDefaultWidth / 2;
constexpr int kProcessImageHeight = kDefaultHeight / 2;
constexpr float kDBRangeMin = -60.0;
constexpr int kThumbnailWidth = 120;
constexpr int kThumbnailHeight = 90;
constexpr int kThumbnailX = 10;
constexpr int kThumbnailY = 10;

template<class PixelType>
struct dlib_pixel_type_to_juce_pixel_format;

template<>
struct dlib_pixel_type_to_juce_pixel_format<dlib::rgb_pixel>
{
    static constexpr juce::Image::PixelFormat pixel_format
    = juce::Image::PixelFormat::RGB;
};

template<>
struct dlib_pixel_type_to_juce_pixel_format<dlib::rgb_alpha_pixel>
{
    static constexpr juce::Image::PixelFormat pixel_format
    = juce::Image::PixelFormat::ARGB;
};

template<class PixelType>
struct JuceImageArray2d
{
    using pixel_type = PixelType;

    JuceImageArray2d(juce::Image img)
    {
        assert(img.isValid());
        assert(img.getPixelData()->pixelFormat
               == dlib_pixel_type_to_juce_pixel_format<PixelType>::pixel_format);
        img_ = img;
    }

    void swap(JuceImageArray2d &rhs)
    {
        std::swap(img_, rhs.img_);
    }

public:
    juce::Image img_;
};

template<class PixelType>
long        num_rows      (const JuceImageArray2d<PixelType>& img)
{
    return img.img_.getHeight();
}

template<class PixelType>
long        num_columns   (const JuceImageArray2d<PixelType>& img)
{
    return img.img_.getWidth();
}

// not supported.
// void        set_image_size(      JuceImageArray2d& img, long rows, long cols);

template<class PixelType>
void*       image_data    (      JuceImageArray2d<PixelType>& img)
{
    auto bmp = juce::Image::BitmapData(img.img_, juce::Image::BitmapData::readWrite);
    return bmp.data;
}

template<class PixelType>
const void* image_data    (const JuceImageArray2d<PixelType>& img)
{
    auto bmp = juce::Image::BitmapData(img.img_, juce::Image::BitmapData::readOnly);
    return bmp.data;
}

template<class PixelType>
long        width_step    (const JuceImageArray2d<PixelType>& img)
{
    auto bmp = juce::Image::BitmapData(img.img_, juce::Image::BitmapData::readOnly);
    return bmp.lineStride;
}

template<class PixelType>
void        swap          (      JuceImageArray2d<PixelType>& a, JuceImageArray2d<PixelType>& b)
{
    a.swap(b);
}

namespace dlib
{
    template <class PixelType>
    struct image_traits<JuceImageArray2d<PixelType>>
    {
        typedef typename JuceImageArray2d<PixelType>::pixel_type pixel_type;
    };
}

template<class T = void>
struct default_div
{
    T operator()(T const &t, double x) const {
        assert(x != 0);
        return t / x;
    }
};

template<>
struct default_div<void>
{
    template<class T>
    T operator()(T const &t, double x) const {
        assert(x != 0);
        return t / x;
    }
};

template<class T, class PlusFunc = std::plus<T>, class DivFunc = default_div<T>>
struct MovingAverageValue
{
    MovingAverageValue()
    :   plus_(PlusFunc{})
    ,   div_(DivFunc{})
    {}

    MovingAverageValue(PlusFunc plus, DivFunc div)
    :   plus_(plus)
    ,   div_(div)
    {}

    int GetNumHistory() const
    {
        return history_.size();
    }

    void SetNumHistory(int n, T const &init = T{})
    {
        history_.resize(n, init);
    }

    template<class U>
    void Push(U &&x)
    {
        history_[write_pos_] = std::move(x);
        write_pos_ = (write_pos_ + 1) % history_.size();
    }

    T GetAverage() const
    {
        T tmp = T{};

        auto const num = history_.size();
        for(int i = 0; i < num; ++i) {
            tmp = plus_(tmp, history_[(i + write_pos_) % num]);
        }

        return div_(tmp, (double)num);
    }

private:
    std::vector<T> history_;
    int write_pos_ = 0;
    PlusFunc plus_;
    DivFunc div_;
};

struct AudioPluginAudioProcessorEditor::Impl
:   public juce::CameraDevice::Listener
,   private juce::Thread
{
    using WF = juce::dsp::WindowingFunction<float>;
    using AP = AudioPluginAudioProcessor;

    static constexpr int kMouthPointsBegin = 48;
    static constexpr int kMouthPointsEnd = 68;
    static constexpr int kNumMouthPoints = kMouthPointsEnd - kMouthPointsBegin;
    static constexpr int kNumHistory = 3;

    Impl(AudioPluginAudioProcessorEditor *owner)
    :   juce::Thread("Image Processing Thread")
    ,   owner_(owner)
    {
        fft_ = std::make_unique<juce::dsp::FFT>(AP::kOrder);
        fft_buffer_.resize(AP::kNumBins * 2);
        wnd_ = std::make_unique<WF>(AP::kNumBins, WF::blackman);

        mtx_call_async_ = std::make_shared<std::mutex>();
        face_data_.SetNumHistory(kNumHistory);
    }

    ~Impl()
    {
        stopImageProcessingThread();
        closeCamera();
    }

    void startImageProcessingThread()
    {
        juce::Thread::startThread();
    }

    void stopImageProcessingThread()
    {
        if(isThreadRunning() == false) { return; }

        auto p = mtx_call_async_;
        {
            assert(p);
            // call_async 処理中は、その終了を待機する。
            std::unique_lock lock(*p);

            mtx_call_async_.reset();
        }

        signalThreadShouldExit();
        juce::Thread::stopThread(5000);
    }

    struct FaceData
    {
        using MouthPointArray = std::array<MovingAverageValue<dlib::dpoint>, kNumMouthPoints>;
        juce::Image image_;             //!< エディター画面で表示するカメラ画像
        MouthPointArray mouth_points_;  //!< エディター画面で表示する口の輪郭
        MovingAverageValue<double> mar_;
        int num_skipped_ = 0;           //!< 顔検出できなかったフレーム数。（検出できた時点で 0 にリセットされる）
        int num_processed_ = 0;

        void SetNumHistory(int n) {
            for(auto &p: mouth_points_) { p.SetNumHistory(n); }
            mar_.SetNumHistory(n);
        }
    };

    FaceData getFaceData() const
    {
        std::unique_lock lock(gui_mtx_);
        return face_data_;
    }

    //! 対数振幅スペクトルを取得
    //! s には、 [0..fs/2] の範囲の最新の対数振幅スペクトルのデータが書き込まれる。
    using SpectrumType = std::array<float, AP::kNumBins/2+1>;

    SpectrumType const & calculateLogAmplitudeSpectrum()
    {
        std::fill(fft_buffer_.begin(), fft_buffer_.end(), 0.0);
        owner_->processorRef.getSampleHistory(fft_buffer_);

        wnd_->multiplyWithWindowingTable(fft_buffer_.data(), AP::kNumBins);

        fft_->performRealOnlyForwardTransform(fft_buffer_.data());

        int const N = AP::kNumBins;
        assert(N == fft_->getSize());

        //! 直流成分とナイキスト周波数成分を 1 / N で正規化する。
        fft_buffer_[0] /= N;
        fft_buffer_[AP::kNumBins] /= N;

        //! それ以外の交流成分は、 1/ N で正規化後、 複素共役の分を考慮して 2 倍する。
        std::transform(fft_buffer_.data() + 1,
                       fft_buffer_.data() + N,
                       fft_buffer_.data() + 1,
                       [&](float f) -> float { return f / (double)N * 2.0; });

        //! 対数振幅スペクトルを計算
        for(int i = 0, end = spectrum_.size(); i < end ; ++i) {
            auto amp_spec = std::sqrt(fft_buffer_[i*2] * fft_buffer_[i*2] +
                                      fft_buffer_[i*2+1] * fft_buffer_[i*2+1]);

            if(amp_spec == 0) {
                spectrum_[i] = kDBRangeMin;
            } else {
                spectrum_[i] = std::max<float>(kDBRangeMin, std::log10(amp_spec) * 20);
            }
        }

        return spectrum_;
    }

    bool openCamera(int index)
    {
        closeCamera();

        // This function crashes on Bitwig Studio 3.1.
        // Is there any workaround?
        cam_.reset(juce::CameraDevice::openDevice(index, kMinWidth, kMinHeight, kMaxWidth, kMaxHeight));
        if(cam_) {
            cam_->addListener(this);

            // このウィンドウは表示しない。
            // （画面が顔の動きと反転しているのでこの画像はそのまま使わない。代わりに画像処理スレッドで処理した画像を反転して表示する）
            viewer_ = cam_->createViewerComponent();
            return true;
        } else {
            return false;
        }
    }

    void closeCamera()
    {
        if(cam_) {
            delete viewer_;
            viewer_ = nullptr;
            cam_->removeListener(this);
            cam_.reset();
        }
    }

    void setImageFlipFlag(bool flip)
    {
        flip_ = flip;
    }

    bool shouldFlipImage() const { return flip_.load(); }

public:
    AudioPluginAudioProcessorEditor *owner_ = nullptr;
    std::unique_ptr<juce::CameraDevice> cam_;
    juce::Component *viewer_ = nullptr;

    juce::ComboBox cmb_camera_list_;
    juce::Slider sl_low_freq_;
    juce::Slider sl_high_freq_;
    juce::Slider sl_qfactor_;

    juce::Label lbl_low_freq_;
    juce::Label lbl_high_freq_;
    juce::Label lbl_qfactor_;

private:
    std::mutex mutable gui_mtx_;    //!< GUI スレッドと画像処理スレッドの排他制御

    std::vector<float> fft_buffer_;
    std::unique_ptr<juce::dsp::FFT> fft_;
    std::unique_ptr<juce::dsp::WindowingFunction<float>> wnd_;

    FaceData face_data_;
    SpectrumType spectrum_;
    std::atomic<bool> flip_ = { true };

    //! callAsync 呼び出し中に確保するべきミューテックス。
    /* AudioPluginAudioProcessorEditor デストラクト時に、
     * 画像処理スレッドを停止したあとで callAsync の処理が呼び出されてアクセス違反が発生することがある。
     * これを回避するため、　callAsync　の中でこのミューテックスを確保することとし、
     * このミューテックスが破棄されているときは何も処理をしないようにする。
     */
    std::shared_ptr<std::mutex> mtx_call_async_;

private:
    void setFaceData(FaceData const &fd)
    {
        std::unique_lock lock(gui_mtx_);
        face_data_ = fd;
    }

    double getMouthAspectRatio(std::vector<dlib::dpoint> const &points)
    {
        auto distance = [](dlib::dpoint pt1, dlib::dpoint pt2) {
            auto const x = pt2.x() - pt1.x();
            auto const y = pt2.y() - pt1.y();
            return std::sqrt(x * x + y * y);
        };

        auto const A = distance(points[2], points[10]);
        auto const B = distance(points[4], points[8]);
        auto const C = distance(points[0], points[6]);

        auto const mar = (A + B) / (2.0 * C);

        return mar;
    }

    struct ImageProcessingContext {
        ImageProcessingContext()
        {
            detector = dlib::get_frontal_face_detector();

            juce::File executable_dir = juce::File::getSpecialLocation(juce::File::currentExecutableFile).getParentDirectory();
            juce::File landmark_file = executable_dir.getChildFile("Data/shape_predictor_68_face_landmarks.dat");
            if(landmark_file.existsAsFile() == false) {
                landmark_file = executable_dir.getParentDirectory().getChildFile("Resources/shape_predictor_68_face_landmarks.dat");
            }

#if JUCE_WINDOWS
            std::ifstream ifs(landmark_file.getFullPathName().toWideCharPointer(), std::ios::binary | std::ios::in);
#else
            std::ifstream ifs(landmark_file.getFullPathName().toUTF8(), std::ios::binary | std::ios::in);
#endif

            dlib::deserialize(predictor, ifs);

            tmp_mouth_points.resize(kNumMouthPoints);
            set_image_size(tmp_array, kProcessImageHeight, kProcessImageWidth);
        }

        dlib::frontal_face_detector detector;
        dlib::shape_predictor predictor;

        std::vector<dlib::dpoint> tmp_mouth_points;
        dlib::array2d<unsigned char> tmp_array;
    };

    ImageProcessingContext context_;

    //! camera から取得した画像を、画像処理スレッドで処理するために保持しておく変数
    juce::Image tmp_image_;

    void ProcessImage(juce::Image img, ImageProcessingContext &ctx)
    {
        static bool write_image = false;
        auto get_desktop_file = [](juce::String const &filename) {
            return juce::File::getSpecialLocation(juce::File::userDesktopDirectory).getChildFile(filename);
        };

        if(write_image) {
            juce::PNGImageFormat fmt;
            juce::File file = get_desktop_file("original_image.png");
            auto fs = file.createOutputStream();
            fmt.writeImageToStream(img, *fs);
            fs.reset();
        }

        // expected image ascpected ratio
        double const ar = (double)kDefaultWidth / kDefaultHeight;

        auto const actual_w = img.getWidth();
        auto const expected_w = (int)std::round(img.getHeight() * ar);

        if(actual_w == expected_w) {
            // do nothing.
        } else if(actual_w > expected_w) {
            // 画面が横長なので左右を捨てる。
            int const nx = actual_w - expected_w;
            auto rect = img.getBounds().reduced(nx / 2.0, 0);
            img = img.getClippedImage(rect);
        } else {
            // 画面が縦長なので上下を捨てる。
             int const ny = img.getHeight() - (img.getWidth() / ar);
             auto rect = img.getBounds().reduced(0, ny / 2.0);
             img = img.getClippedImage(rect);
        }

        img = img.rescaled(kProcessImageWidth, kProcessImageHeight);
        img.desaturate();

        if(write_image) {
            juce::PNGImageFormat fmt;
            juce::File file = get_desktop_file("juce_image.png");
            auto fs = file.createOutputStream();
            fmt.writeImageToStream(img, *fs);
            fs.reset();
        }

        if(img.getPixelData()->pixelFormat == juce::Image::PixelFormat::RGB) {
            dlib::assign_image(ctx.tmp_array, JuceImageArray2d<dlib::rgb_pixel>(img));
        } else {
            dlib::assign_image(ctx.tmp_array, JuceImageArray2d<dlib::rgb_alpha_pixel>(img));
        }

#ifdef DLIB_PNG_SUPPORT
        if(write_image) {
            dlib::save_png(ctx.tmp_array, get_desktop_file("dlib_image.png").getFullPathName().toStdString());
        }
#endif

        auto rects = ctx.detector(ctx.tmp_array, 0);

        FaceData tmp_face_data = getFaceData();

        if(rects.empty()) {
            tmp_face_data.num_skipped_ += 1;
            if(kShowFace) {
                tmp_face_data.image_ = std::move(img);
            }

            setFaceData(tmp_face_data);
        } else {
            dlib::full_object_detection shape = ctx.predictor(ctx.tmp_array, rects[0]);

            for(int i = 0, end = kNumMouthPoints; i < end; ++i) {
                ctx.tmp_mouth_points[i] = shape.part(i + kMouthPointsBegin);
            }

            double const tmp_mar = getMouthAspectRatio(ctx.tmp_mouth_points);

            for(int i = 0, end = kNumMouthPoints; i < end; ++i) {
                tmp_face_data.mouth_points_[i].Push(ctx.tmp_mouth_points[i]);
            }

            tmp_face_data.mar_.Push(tmp_mar);
            tmp_face_data.num_processed_ += 1;
            tmp_face_data.num_skipped_ = 0;

            if(kShowFace) {
                tmp_face_data.image_ = std::move(img);
            }

            setFaceData(tmp_face_data);
        }
    }

    void run() override
    {
        int continue_count_ = 0;

        for( ; threadShouldExit() == false; ) {

            juce::Image tmp_image = std::move(tmp_image_);

            if(tmp_image.isValid() == false) {
                continue_count_ += 1;
                if(continue_count_ < 10) {
                    // do nothing.
                } else if(continue_count_ < 30) {
                    std::this_thread::yield();
                } else {
                    std::this_thread::sleep_for(std::chrono::milliseconds(10));
                }

                continue;
            }

            continue_count_ = 0;

            ProcessImage(tmp_image, context_);

            auto mm = juce::MessageManager::getInstance();
            mm->callAsync([this, w = std::weak_ptr<std::mutex>(mtx_call_async_)] {

                for( ; ; ) {
                    auto const pmtx = w.lock();
                    if(!pmtx) {
                        // callAsync 呼び出し中に確保するべき mutex が破棄されている。
                        // => 何も処理をせずにリターンする。
                        return;
                    }

                    std::unique_lock lock(*pmtx, std::try_to_lock);
                    if(lock) {
                        owner_->OnUpdateMouthAspectRatio();
                        break;
                    }

                    std::this_thread::yield();
                }
            });
        }
    }

    void imageReceived(juce::Image const &image) override
    {
#if WAHWTH_USE_IMAGE_PROCESSING_THREAD
        tmp_image_ = image;
#else
        ProcessImage(image, context_);

        auto mm = juce::MessageManager::getInstance();
        mm->callAsync([this] {
            owner_->OnUpdateMouthAspectRatio();
        });
#endif
    }
};

//==============================================================================
AudioPluginAudioProcessorEditor::AudioPluginAudioProcessorEditor (AudioPluginAudioProcessor& p)
:   AudioProcessorEditor(&p)
,   processorRef(p)
,   pimpl_(std::make_unique<Impl>(this))
{
#if JUCE_WINDOWS
    String typeFaceName = "Meiryo UI";
    Desktop::getInstance().getDefaultLookAndFeel().setDefaultSansSerifTypefaceName(typeFaceName);
#elif JUCE_MAC
    String typeFaceName = "Arial Unicode MS";
    Desktop::getInstance().getDefaultLookAndFeel().setDefaultSansSerifTypefaceName(typeFaceName);
#endif

    juce::ignoreUnused (processorRef);

    setResizable(true, true);
    setResizeLimits(kMinWidth, kMinHeight, kMaxWidth, kMaxHeight);

    // Make sure that before the constructor has finished, you've set the
    // editor's size to whatever you need it to be.
    setSize (kDefaultWidth, kDefaultHeight);

    addAndMakeVisible(pimpl_->cmb_camera_list_);
    addAndMakeVisible(pimpl_->sl_qfactor_);
    addAndMakeVisible(pimpl_->sl_low_freq_);
    addAndMakeVisible(pimpl_->sl_high_freq_);
    addAndMakeVisible(pimpl_->lbl_low_freq_);
    addAndMakeVisible(pimpl_->lbl_high_freq_);
    addAndMakeVisible(pimpl_->lbl_qfactor_);

    pimpl_->cmb_camera_list_.addItemList(juce::CameraDevice::getAvailableDevices(), 1);
    pimpl_->sl_qfactor_.setSliderStyle(juce::Slider::LinearBar);
    pimpl_->sl_low_freq_.setSliderStyle(juce::Slider::LinearBar);
    pimpl_->sl_high_freq_.setSliderStyle(juce::Slider::LinearBar);

    pimpl_->sl_low_freq_.setTextValueSuffix(" Hz");
    pimpl_->sl_high_freq_.setTextValueSuffix(" Hz");

    auto make_value_from_text_function = [](juce::Slider *slider) {
        return [slider](juce::String const &str) -> double {
            auto range = slider->getRange();
            try {
                return std::clamp<double>(std::stof(str.toStdString()), range.getStart(), range.getEnd());
            } catch(...) {
                return range.getStart();
            }
        };
    };

    auto make_text_from_value_function = [](juce::Slider *slider) {
        return [slider](double value) -> juce::String {
            char buf[32] = {};
            sprintf(buf, "%0.2f", value);
            return buf;
        };
    };

    pimpl_->sl_low_freq_.valueFromTextFunction = make_value_from_text_function(&pimpl_->sl_low_freq_);
    pimpl_->sl_low_freq_.textFromValueFunction = make_text_from_value_function(&pimpl_->sl_low_freq_);
    pimpl_->sl_high_freq_.valueFromTextFunction = make_value_from_text_function(&pimpl_->sl_high_freq_);
    pimpl_->sl_high_freq_.textFromValueFunction = make_text_from_value_function(&pimpl_->sl_high_freq_);
    pimpl_->sl_qfactor_.valueFromTextFunction = make_value_from_text_function(&pimpl_->sl_qfactor_);
    pimpl_->sl_qfactor_.textFromValueFunction = make_text_from_value_function(&pimpl_->sl_qfactor_);

    using AP = AudioPluginAudioProcessor;
    pimpl_->sl_low_freq_.setRange(AP::kLowFreqMin, AP::kLowFreqMax);
    pimpl_->sl_high_freq_.setRange(AP::kHighFreqMin, AP::kHighFreqMax);
    pimpl_->sl_qfactor_.setRange(AP::kQFactorMin, AP::kQFactorMax);

    pimpl_->sl_low_freq_.setValue(*processorRef.low_freq_, juce::dontSendNotification);
    pimpl_->sl_high_freq_.setValue(*processorRef.high_freq_, juce::dontSendNotification);
    pimpl_->sl_qfactor_.setValue(*processorRef.qfactor_, juce::dontSendNotification);

    pimpl_->sl_low_freq_.setSkewFactor(0.6);
    pimpl_->sl_high_freq_.setSkewFactor(0.6);
    pimpl_->sl_qfactor_.setSkewFactor(1.0);

    pimpl_->cmb_camera_list_.addListener(this);
    pimpl_->sl_low_freq_.addListener(this);
    pimpl_->sl_high_freq_.addListener(this);
    pimpl_->sl_qfactor_.addListener(this);

    pimpl_->lbl_low_freq_.setText ("Freq Lo", juce::NotificationType::dontSendNotification);
    pimpl_->lbl_high_freq_.setText("Freq Hi", juce::NotificationType::dontSendNotification);
    pimpl_->lbl_qfactor_.setText  ("Resonance", juce::NotificationType::dontSendNotification);

    auto ed = processorRef.editor_data_.load();

    auto const num_cameras = pimpl_->cmb_camera_list_.getNumItems();
    if(num_cameras != 0) {
        // 保存されたカメラ番号
        auto const try_first = std::clamp(ed.camera_index_,
                                          AP::kCameraIndexMin,
                                          num_cameras);

        // カメラのオープンを試す順序を用意する。
        // 保存されたカメラ番号を最初に試す。
        // 一番最初にオープンできたカメラを使用する。
        std::vector<int> camera_open_order;
        camera_open_order.push_back(try_first);
        for(int i = 0, end = num_cameras; i < end; ++i) {
            if(i != try_first) { camera_open_order.push_back(i); }
        }

        for(auto index : camera_open_order) {
            if(pimpl_->openCamera(index)) {
                pimpl_->cmb_camera_list_.setSelectedItemIndex(index,
                                                              juce::NotificationType::dontSendNotification
                                                              );
                break;
            }
        }
    }

    pimpl_->setImageFlipFlag(ed.flip_);

    processorRef.SetEditorDataUpdateCallback([this]() {
        auto ed = processorRef.editor_data_.load();

        auto const index = std::clamp<int>(ed.camera_index_,
                                           AP::kCameraIndexMin,
                                           pimpl_->cmb_camera_list_.getNumItems());

        pimpl_->cmb_camera_list_.setSelectedItemIndex(index);
        pimpl_->setImageFlipFlag(ed.flip_);
    });

    getConstrainer()->setFixedAspectRatio((double)kMaxWidth / kMaxHeight);
    setBounds(0, 0, kDefaultWidth, kDefaultHeight);
    resized();

#if WAHWTH_USE_IMAGE_PROCESSING_THREAD
    pimpl_->startImageProcessingThread();
#endif
}

AudioPluginAudioProcessorEditor::~AudioPluginAudioProcessorEditor()
{
    processorRef.SetEditorDataUpdateCallback(nullptr);

    pimpl_->closeCamera();

#if WAHWTH_USE_IMAGE_PROCESSING_THREAD
    pimpl_->stopImageProcessingThread();
#endif

    pimpl_.reset();
}

void AudioPluginAudioProcessorEditor::OnUpdateMouthAspectRatio()
{
    auto const face_data = pimpl_->getFaceData();

    double const kMarMin = 0.4;
    double const kMarMax = 0.95;

    double const mar = std::clamp(face_data.mar_.GetAverage(), kMarMin, kMarMax);
    double const freq = (mar - kMarMin) / (kMarMax - kMarMin);

    auto *param = processorRef.freq_;
    param->setValueNotifyingHost(freq);

    repaint();
}

void AudioPluginAudioProcessorEditor::sliderValueChanged(juce::Slider *s)
{
    if(s == &pimpl_->sl_qfactor_) {
        auto *param = processorRef.qfactor_;
        param->setValueNotifyingHost(param->convertTo0to1(s->getValue()));
    } else if(s == &pimpl_->sl_low_freq_) {
        auto *param = processorRef.low_freq_;
        param->setValueNotifyingHost(param->convertTo0to1(s->getValue()));
    } else if(s == &pimpl_->sl_high_freq_) {
        auto *param = processorRef.high_freq_;
        param->setValueNotifyingHost(param->convertTo0to1(s->getValue()));
    } else {
        assert(false);
    }
}

void AudioPluginAudioProcessorEditor::comboBoxChanged(juce::ComboBox *c)
{
    if(c == &pimpl_->cmb_camera_list_) {
        auto const item_index = c->getSelectedItemIndex();
        if(pimpl_->openCamera(item_index)) {
            auto tmp = processorRef.editor_data_.load();
            tmp.camera_index_ = item_index;
            processorRef.editor_data_ = tmp;
        } else {
            c->setSelectedId(0, juce::NotificationType::dontSendNotification);
        }
    } else {
        assert(false);
    }
}

//==============================================================================
void AudioPluginAudioProcessorEditor::paint (juce::Graphics& g)
{
    // (Our component is opaque, so we must completely fill the background with a solid colour)
    g.fillAll (getLookAndFeel().findColour (juce::ResizableWindow::backgroundColourId));

    auto const face_data = pimpl_->getFaceData();

    if(kShowFace && face_data.image_.isValid()) {
        double scale = face_data.image_.getWidth() / (double)kThumbnailWidth;
        juce::AffineTransform t;
        if(pimpl_->shouldFlipImage()) {
            t = t.scaled(1.0 / scale);
            t = t.scaled(-1, 1.0); // カメラから取得した画像を、顔の動きに合うように左右反転して表示する。
            t = t.translated(kThumbnailWidth + kThumbnailX, kThumbnailY);
        } else {
            t = t.scaled(1.0 / scale);
            t = t.translated(kThumbnailX, kThumbnailY);
        }
        g.drawImageTransformed(face_data.image_, t);
     }

    g.setFont (15.0f);

    auto &points = face_data.mouth_points_;

    auto const col_bright = juce::Colour(0.3f, 0.7f, 0.9f, 0.75f);
    auto const col_dark = juce::Colour(0.3f, 0.8f, 0.2f, 0.75f);

    g.setColour(col_bright);

    // 口の形状を描画する
    if(face_data.num_skipped_ <= 2) {
        int left = 10000, top = 10000, right = 0, bottom = 0;

        for(auto const &pts: points) {
            auto avg = pts.GetAverage();
            left = std::min<int>(left, avg.x());
            right = std::max<int>(right, avg.x());
            top = std::min<int>(top, avg.y());
            bottom = std::max<int>(bottom, avg.y());
        }

        juce::Rectangle<int> b_mouth {
            left,
            top,
            std::max<int>(1, (right - left)), // 0除算を回避
            std::max<int>(1, (bottom - top))  // 0除算を回避
        };

        juce::Rectangle<int> b_area = getBounds();

        b_area = b_area.reduced((int)(b_area.getWidth() * 0.2), (int)(b_area.getHeight() * 0.2));

        auto const shift_mouth = b_mouth.getCentre();
        auto const shift_area = b_area.getCentre();

        double const scale_x = b_area.getWidth() / (double)b_mouth.getWidth();
        double const scale_y = b_area.getHeight() / (double)b_mouth.getHeight();
        double const scale = std::min<double>(scale_x, scale_y);

        b_mouth = b_mouth.expanded(b_mouth.getWidth() * (scale - 1.0) / 2.0,
                                   b_mouth.getHeight() * (scale - 1.0) / 2.0);

        auto drawLine = [&](int a, int b) {
            auto avg_a = points[a].GetAverage();
            auto avg_b = points[b].GetAverage();

            auto ax = avg_a.x();
            auto ay = avg_a.y();
            auto bx = avg_b.x();
            auto by = avg_b.y();

            if(pimpl_->shouldFlipImage()) {
                ax = (ax - b_mouth.getCentre().x) * -1 + b_mouth.getCentre().x;
                bx = (bx - b_mouth.getCentre().x) * -1 + b_mouth.getCentre().x;
            }

            g.drawLine((ax - shift_mouth.x) * scale + shift_area.x,
                       (ay - shift_mouth.y) * scale + shift_area.y,
                       (bx - shift_mouth.x) * scale + shift_area.x,
                       (by - shift_mouth.y) * scale + shift_area.y
                       );
        };

        for(int i = 0; i < 12; ++i) {
            drawLine(i, (i+1)%12);
        }

        for(int i = 0; i < 8; ++i) {
            drawLine(i+12, ((i+1)%8)+12);
        }
    }

    // スペクトルを描画する
    auto const &spectrum = pimpl_->calculateLogAmplitudeSpectrum();
    auto b = getBounds().reduced(10);

    b.removeFromTop(90);

    g.setColour(col_bright);

    static double kLogEmphasisCoeff = 200; // この値が大きいほうが低音域が広く表示される。 [9 .. inf]
    float last_x = -1;
    for(int i = 0, end = spectrum.size()-1; i < end; ++i) {
        auto &s1 = spectrum[i];
        auto &s2 = spectrum[i+1];

        float x1 = log10((i / (double)end) * kLogEmphasisCoeff + 1.0) / log10(kLogEmphasisCoeff + 1.0) * b.getWidth() + b.getX();
        float x2 = log10(((i+1) / (double)end) * kLogEmphasisCoeff + 1.0) / log10(kLogEmphasisCoeff + 1.0) * b.getWidth() + b.getX();
        float y1 = (s1 / (double)kDBRangeMin) * b.getHeight() + b.getY();
        float y2 = (s2 / (double)kDBRangeMin) * b.getHeight() + b.getY();

        assert(last_x < x1);
        assert(x1 < x2);
        g.drawLine(x1, y1, x2, y2, 2.0);
    }

    char const nl = '\n';
    std::stringstream ss;

    ss
    << "Mouth Aspect Ratio : " << face_data.mar_.GetAverage();

    g.drawFittedText (ss.str(), getLocalBounds().removeFromTop(100), juce::Justification::centred, 1);
}

void AudioPluginAudioProcessorEditor::resized()
{
    // This is generally where you'll want to lay out the positions of any
    // subcomponents in your editor..

    auto b = getBounds().reduced(10);
    auto right = b.removeFromRight(b.getWidth() / 3);
    pimpl_->cmb_camera_list_.setBounds(right.removeFromTop(20));
    right.removeFromTop(2);

    auto align_label_and_control = [](juce::Component &label,
                                      juce::Component &control,
                                      juce::Rectangle<int> b) {
        label.setBounds(b.removeFromLeft((int)(b.getWidth() * 0.33)));
        control.setBounds(b);
    };

    align_label_and_control(pimpl_->lbl_low_freq_, pimpl_->sl_low_freq_, right.removeFromTop(20));
    right.removeFromTop(2);

    align_label_and_control(pimpl_->lbl_high_freq_, pimpl_->sl_high_freq_, right.removeFromTop(20));
    right.removeFromTop(2);

    align_label_and_control(pimpl_->lbl_qfactor_, pimpl_->sl_qfactor_, right.removeFromTop(20));
}

void AudioPluginAudioProcessorEditor::mouseUp(juce::MouseEvent const&e)
{
    auto const p = e.getPosition();

    // サムネイル領域で右クリックした場合は画像の反転メニューを表示する。
    if(e.mods == juce::ModifierKeys::rightButtonModifier &&
       kThumbnailX <= p.x && p.x < (kThumbnailX + kThumbnailWidth) &&
       kThumbnailY <= p.y && p.y < (kThumbnailY + kThumbnailHeight))
    {
        juce::PopupMenu m;
        m.addItem("Invert", [this]() {
            auto const new_flag = pimpl_->shouldFlipImage() == false;
            pimpl_->setImageFlipFlag(new_flag);
            auto tmp = processorRef.editor_data_.load();
            tmp.flip_ = new_flag;
            processorRef.editor_data_ = tmp;
        });

        m.show();
        return;
    }

    // それ以外の場所を右クリックした場合は About ダイアログを表示する。
    if(e.mods == juce::ModifierKeys::rightButtonModifier) {
        juce::PopupMenu menu;
        menu.addItem(1, "About...");

        int result = menu.show();
        if(result == 1) {
            showModalAboutDialog(this);
        }
    }
}

max_stars_repo_name: hotwatermorning/wahwth

max_issues_repo_name: hotwatermorning/wahwth

max_forks_repo_name: hotwatermorning/wahwth


text:
#include <cradle/gui/background_job_request.hpp>

#include <cradle/background/internals.hpp>
#include <cradle/gui/internals.hpp>
#include <cradle/gui/app/internals.hpp>

#include <boost/thread/thread.hpp>

namespace cradle {

bool update_general_background_job(
    gui_context& ctx, background_job_data& data,
    id_interface const& id,
    boost::function<background_job_interface*()> const& create_background_job)
{
    assert(is_refresh_pass(ctx));

    bool changed = false;

    // If the request is gettable, but the pointer isn't initialized or
    // doesn't have the same ID, reset it to the new request.
    if (!data.ptr.is_initialized() ||
        data.ptr.key() != id)
    {
        auto const& bg = ctx.gui_system->bg;
        data.ptr.reset(*bg, id);
        add_untyped_background_job(data.ptr, *bg,
            // I think this is a legacy interface, and we don't know what the
            // job is going to do, so mark it as a write job to be safe.
            background_job_queue_type::WEB_WRITE,
            create_background_job());
        inc_version(data.abbreviated_identity);
        changed = true;
    }

    // If we already have the result, we're done.
    if (data.ptr.is_ready())
        return changed;

    // Otherwise, update to bring in changes from the background.
    data.ptr.update();

    // Check again to see if that made the pointer ready.
    if (data.ptr.is_ready())
        changed = true;

    request_refresh(ctx, 1);

    return changed;
}

}

max_stars_repo_name: dotdecimal/open-cradle

max_issues_repo_name: dotdecimal/open-cradle

max_forks_repo_name: dotdecimal/open-cradle


text:
//  boost/chrono/stopwatches/reporters/laps_stopwatch_default_formatter.hpp
//  Copyright 2011 Vicente J. Botet Escriba
//  Copyright (c) Microsoft Corporation 2014
//  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://www.boost.org/libs/chrono/stopwatches for documentation.

#ifndef BOOST_CHRONO_STOPWATCHES_REPORTERS_LAPS_STOPWATCH_DEFAULT_FORMATTER_HPP
#define BOOST_CHRONO_STOPWATCHES_REPORTERS_LAPS_STOPWATCH_DEFAULT_FORMATTER_HPP

#include <boost/chrono/stopwatches/reporters/stopwatch_reporter_default_formatter.hpp>
#include <boost/chrono/stopwatches/formatters/accumulator_set_formatter.hpp>
#include <boost/chrono/stopwatches/collectors/laps_accumulator_set.hpp>
#include <boost/chrono/stopwatches/laps_stopwatch.hpp>
#if ! BOOST_OS_WINDOWS || BOOST_PLAT_WINDOWS_DESKTOP
#include <boost/chrono/thread_clock.hpp>
#endif

namespace boost
{
  namespace chrono
  {

    template <typename CharT, typename Clock, typename Features, typename Weight>
    struct basic_stopwatch_reporter_default_formatter<CharT, laps_stopwatch<Clock,laps_accumulator_set<typename Clock::duration,Features, Weight> > >
    {
      typedef basic_accumulator_set_formatter<milli,CharT> type;
    };

  } // namespace chrono
} // namespace boost

#endif

max_stars_repo_name: EdwardPrentice/wrongo

max_issues_repo_name: wujf/mongo

max_forks_repo_name: wujf/mongo


text:
// Copyright 2008-2022 Emil Dotchevski and Reverge Studios, Inc.

// 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)

#ifndef BOOST_QVM_VEC_TRAITS_HPP_INCLUDED
#define BOOST_QVM_VEC_TRAITS_HPP_INCLUDED

#include <boost/qvm/is_scalar.hpp>
#include <boost/qvm/enable_if.hpp>
#include <boost/qvm/config.hpp>

namespace boost { namespace qvm {

template <class V>
struct
vec_traits
    {
    static int const dim=0;
    typedef void scalar_type;
    };

template <class T>
struct
is_vec
    {
    static bool const value = is_scalar<typename vec_traits<T>::scalar_type>::value && vec_traits<T>::dim>0;
    };

namespace
qvm_detail
    {
    template <class T, T>
    struct
    vtr_dispatch_yes
        {
        char x, y;
        };
    }

template <class T>
class
vec_write_element_ref
    {
    template <class U>
    static qvm_detail::vtr_dispatch_yes<typename vec_traits<U>::scalar_type & (*)( U & ), &vec_traits<U>::template write_element<0> > check(int);

    template <class>
    static char check(long);

    public:

    static bool const value = sizeof(check<T>(0)) > 1;
    };

template <int I, class V>
BOOST_QVM_CONSTEXPR BOOST_QVM_INLINE_CRITICAL
typename enable_if_c<
    vec_write_element_ref<V>::value,
    void>::type
write_vec_element( V & v, typename vec_traits<V>::scalar_type s )
    {
    vec_traits<V>::template write_element<I>(v) = s;
    }

template <int I, class V>
BOOST_QVM_CONSTEXPR BOOST_QVM_INLINE_CRITICAL
typename enable_if_c<
    !vec_write_element_ref<V>::value,
    void>::type
write_vec_element( V & v, typename vec_traits<V>::scalar_type s )
    {
    vec_traits<V>::template write_element<I>(v, s);
    }

template <class V>
BOOST_QVM_CONSTEXPR BOOST_QVM_INLINE_CRITICAL
typename enable_if_c<
    vec_write_element_ref<V>::value,
    void>::type
write_vec_element_idx( int i, V & v, typename vec_traits<V>::scalar_type s )
    {
    vec_traits<V>::write_element_idx(i, v) = s;
    }

template <class V>
BOOST_QVM_CONSTEXPR BOOST_QVM_INLINE_CRITICAL
typename enable_if_c<
    !vec_write_element_ref<V>::value,
    void>::type
write_vec_element_idx( int i, V & v, typename vec_traits<V>::scalar_type s )
    {
    vec_traits<V>::write_element_idx(i, v, s);
    }

} }

#endif

max_stars_repo_name: boostorg/boost-qvm

max_issues_repo_name: boostorg/boost-qvm

max_forks_repo_name: boostorg/boost-qvm


text:
// Copyright (C) 2010-2015 Joshua Boyce
// See the file COPYING for copying permission.

#include <hadesmem/thread_list.hpp>
#include <hadesmem/thread_list.hpp>

#include <utility>

#include <hadesmem/detail/warning_disable_prefix.hpp>
#include <boost/detail/lightweight_test.hpp>
#include <hadesmem/detail/warning_disable_suffix.hpp>

#include <hadesmem/config.hpp>
#include <hadesmem/detail/smart_handle.hpp>
#include <hadesmem/process_entry.hpp>
#include <hadesmem/thread.hpp>

DWORD WINAPI SleepWrapper(LPVOID /*param*/)
{
  ::Sleep(INFINITE);
  return 0UL;
}

void TestThreadList()
{
  using ThreadListIterCat =
    std::iterator_traits<hadesmem::ThreadList::iterator>::iterator_category;
  HADESMEM_DETAIL_STATIC_ASSERT(
    std::is_base_of<std::input_iterator_tag, ThreadListIterCat>::value);
  using ThreadListConstIterCat = std::iterator_traits<
    hadesmem::ThreadList::const_iterator>::iterator_category;
  HADESMEM_DETAIL_STATIC_ASSERT(
    std::is_base_of<std::input_iterator_tag, ThreadListConstIterCat>::value);

  {
    hadesmem::detail::SmartHandle wait_thread_1(
      ::CreateThread(nullptr, 0, &SleepWrapper, nullptr, 0, nullptr));
    hadesmem::detail::SmartHandle wait_thread_2(
      ::CreateThread(nullptr, 0, &SleepWrapper, nullptr, 0, nullptr));
  }

  hadesmem::ThreadList const thread_list_1(::GetCurrentProcessId());
  hadesmem::ThreadList thread_list_2(thread_list_1);
  hadesmem::ThreadList thread_list_3(std::move(thread_list_2));
  thread_list_2 = std::move(thread_list_3);
  BOOST_TEST(std::begin(thread_list_2) != std::end(thread_list_2));

  auto iter = std::begin(thread_list_1);
  BOOST_TEST(iter != std::end(thread_list_1));
  BOOST_TEST(++iter != std::end(thread_list_1));
  BOOST_TEST_NE(iter->GetId(), 0U);
  DWORD const second_id = iter->GetId();
  BOOST_TEST(++iter != std::end(thread_list_1));
  DWORD const third_id = iter->GetId();
  BOOST_TEST_NE(second_id, third_id);
}

void TestThreadListAlgorithm()
{
  hadesmem::ThreadList const thread_list_1(::GetCurrentProcessId());

  for (auto const& entry : thread_list_1)
  {
    BOOST_TEST_EQ(entry.GetUsage(), 0UL);
    BOOST_TEST_NE(entry.GetId(), 0UL);
    BOOST_TEST_EQ(entry.GetOwnerId(), ::GetCurrentProcessId());
    BOOST_TEST(entry.GetBasePriority() >= 0L);
    BOOST_TEST(entry.GetBasePriority() <= 31L);
    BOOST_TEST_EQ(entry.GetDeltaPriority(), 0L);
    BOOST_TEST_EQ(entry.GetFlags(), 0UL);
  }

  auto const this_iter = std::find_if(std::begin(thread_list_1),
                                      std::end(thread_list_1),
                                      [](hadesmem::ThreadEntry const& thread)
                                      {
    return thread.GetId() == ::GetCurrentThreadId();
  });
  BOOST_TEST(this_iter != std::end(thread_list_1));
}

int main()
{
  TestThreadList();
  TestThreadListAlgorithm();
  return boost::report_errors();
}


text:
#ifndef BOOST_THREAD_TSS_HPP
#define BOOST_THREAD_TSS_HPP
// 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)
// (C) Copyright 2007-8 Anthony Williams

#include <boost/thread/detail/config.hpp>
#include <boost/shared_ptr.hpp>
#include <boost/thread/detail/thread_heap_alloc.hpp>

#include <boost/config/abi_prefix.hpp>

namespace mars_boost {} namespace boost = mars_boost; namespace mars_boost
{
    namespace detail
    {
        struct tss_cleanup_function
        {
            virtual ~tss_cleanup_function()
            {}

            virtual void operator()(void* data)=0;
        };

        BOOST_THREAD_DECL void set_tss_data(void const* key,mars_boost::shared_ptr<tss_cleanup_function> func,void* tss_data,bool cleanup_existing);
        BOOST_THREAD_DECL void* get_tss_data(void const* key);
    }

    template <typename T>
    class thread_specific_ptr
    {
    private:
        thread_specific_ptr(thread_specific_ptr&);
        thread_specific_ptr& operator=(thread_specific_ptr&);

        struct delete_data:
            detail::tss_cleanup_function
        {
            void operator()(void* data)
            {
                delete static_cast<T*>(data);
            }
        };

        struct run_custom_cleanup_function:
            detail::tss_cleanup_function
        {
            void (*cleanup_function)(T*);

            explicit run_custom_cleanup_function(void (*cleanup_function_)(T*)):
                cleanup_function(cleanup_function_)
            {}

            void operator()(void* data)
            {
                cleanup_function(static_cast<T*>(data));
            }
        };

        mars_boost::shared_ptr<detail::tss_cleanup_function> cleanup;

    public:
        typedef T element_type;

        thread_specific_ptr():
            cleanup(detail::heap_new<delete_data>(),detail::do_heap_delete<delete_data>())
        {}
        explicit thread_specific_ptr(void (*func_)(T*))
        {
            if(func_)
            {
                cleanup.reset(detail::heap_new<run_custom_cleanup_function>(func_),detail::do_heap_delete<run_custom_cleanup_function>());
            }
        }
        ~thread_specific_ptr()
        {
            detail::set_tss_data(this,mars_boost::shared_ptr<detail::tss_cleanup_function>(),0,true);
        }

        T* get() const
        {
            return static_cast<T*>(detail::get_tss_data(this));
        }
        T* operator->() const
        {
            return get();
        }
        typename mars_boost::detail::sp_dereference< T >::type operator*() const
        {
            return *get();
        }
        T* release()
        {
            T* const temp=get();
            detail::set_tss_data(this,mars_boost::shared_ptr<detail::tss_cleanup_function>(),0,false);
            return temp;
        }
        void reset(T* new_value=0)
        {
            T* const current_value=get();
            if(current_value!=new_value)
            {
                detail::set_tss_data(this,cleanup,new_value,true);
            }
        }
    };
}

#include <boost/config/abi_suffix.hpp>

#endif

max_stars_repo_name: jonetomtom/mars

max_issues_repo_name: jonetomtom/mars

max_forks_repo_name: pengjinning/mars


text:
//  Copyright John Maddock 2006.
//  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)

#define BOOST_ENABLE_ASSERT_HANDLER
#define BOOST_MATH_MAX_SERIES_ITERATION_POLICY INT_MAX
// for consistent behaviour across compilers/platforms:
#define BOOST_MATH_PROMOTE_DOUBLE_POLICY false
// overflow to infinity is OK, we treat these as zero error as long as the sign is correct!
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error

#include <iostream>
#include <ctime>
#include <boost/multiprecision/mpfr.hpp>
#include <boost/math/special_functions/hypergeometric_1F1.hpp>
#include <boost/math/special_functions/hypergeometric_pFq.hpp>
#include <boost/math/special_functions/relative_difference.hpp>

#include <boost/random.hpp>
#include <set>
#include <fstream>
#include <boost/iostreams/tee.hpp>
#include <boost/iostreams/stream.hpp>

typedef double test_type;
using boost::multiprecision::mpfr_float;

namespace boost {
   //
   // We convert assertions into exceptions, so we can log them and continue:
   //
   void assertion_failed(char const * expr, char const *, char const * file, long line)
   {
      std::ostringstream oss;
      oss << file << ":" << line << " Assertion failed: " << expr;
      throw std::runtime_error(oss.str());
   }

}

void print_value(double x, std::ostream& os = std::cout)
{
   int e;
   double m = std::frexp(x, &e);
   m = std::ldexp(m, 54);
   e -= 54;
   boost::int64_t val = (boost::int64_t)m;
   BOOST_ASSERT(std::ldexp((double)val, e) == x);
   os << "std::ldexp((double)" << val << ", " << e << ")";
}

void print_row(double a, double b, double z, mpfr_float result, std::ostream& os = std::cout)
{
   os << "  {{ ";
   print_value(a, os);
   os << ", ";
   print_value(b, os);
   os << ", ";
   print_value(z, os);
   os << ", SC_(" << std::setprecision(45) << result << ") }}" << std::endl;
}

struct error_data
{
   error_data(double a, double b, double z, boost::intmax_t e)
      : a(a), b(b), z(z), error(e) {}
   double a, b, z;
   boost::intmax_t error;
   bool operator<(const error_data& other)const
   {
      return error < other.error;
   }
};

int main()
{
   try {
      test_type max_a, max_b, max_z, min_a, min_b, min_z;

      unsigned number_of_samples;

      std::ofstream log_stream, incalculable_stream, unevaluated_stream, bins_stream;
      std::string basename;

      std::cout << "Enter range for a: ";
      std::cin >> min_a >> max_a;
      std::cout << "Enter range for b: ";
      std::cin >> min_b >> max_b;
      std::cout << "Enter range for z: ";
      std::cin >> min_z >> max_z;
      std::cout << "Enter number of samples: ";
      std::cin >> number_of_samples;
      std::cout << "Enter basename for log files: ";
      std::cin >> basename;

      typedef boost::iostreams::tee_device<std::ostream, std::ostream> tee_sink;
      typedef boost::iostreams::stream<tee_sink> tee_stream;

      log_stream.open((basename + ".log").c_str());
      tee_stream tee_log(tee_sink(std::cout, log_stream));
      incalculable_stream.open((basename + "_incalculable.log").c_str());
      unevaluated_stream.open((basename + "_unevaluated.log").c_str());
      bins_stream.open((basename + "_bins.csv").c_str());
      tee_stream tee_bins(tee_sink(std::cout, bins_stream));

      boost::random::mt19937 gen(std::time(0));
      boost::random::uniform_real_distribution<test_type> a_dist(min_a, max_a);
      boost::random::uniform_real_distribution<test_type> b_dist(min_b, max_b);
      boost::random::uniform_real_distribution<test_type> z_dist(min_z, max_z);

      std::multiset<error_data> errors;
      std::map<std::pair<int, int>, int> bins;

      unsigned incalculable = 0;
      unsigned evaluation_errors = 0;
      test_type max_error = 0;

      do
      {
         test_type a = a_dist(gen);
         test_type b = b_dist(gen);
         test_type z = z_dist(gen);
         test_type found, expected;
         mpfr_float mp_expected;

         try {
            mp_expected = boost::math::hypergeometric_pFq_precision({ mpfr_float(a) }, { mpfr_float(b) }, mpfr_float(z), 25, 200.0);
            expected = (test_type)mp_expected;
         }
         catch (const std::exception&)
         {
            // Unable to compute reference value:
            ++incalculable;
            tee_log << "Unable to compute reference value in reasonable time: " << std::endl;
            print_row(a, b, z, mpfr_float(0), tee_log);
            incalculable_stream << std::setprecision(6) << std::scientific << a << "," << b << "," << z << "\n";
            continue;
         }
         try
         {
            found = boost::math::hypergeometric_1F1(a, b, z);
         }
         catch (const std::exception&)
         {
            ++evaluation_errors;
            --number_of_samples;
            log_stream << "Unexpected exception calculating value: " << std::endl;
            print_row(a, b, z, mp_expected, log_stream);
            unevaluated_stream << std::setprecision(6) << std::scientific << a << "," << b << "," << z << "\n";
            continue;
         }
         test_type err = boost::math::epsilon_difference(found, expected);
         if (err > max_error)
         {
            tee_log << "New maximum error is: " << err << std::endl;
            print_row(a, b, z, mp_expected, tee_log);
            max_error = err;
         }
         try {
            errors.insert(error_data(a, b, z, boost::math::lltrunc(err)));
         }
         catch (...)
         {
            errors.insert(error_data(a, b, z, INT_MAX));
         }
         --number_of_samples;
         if (number_of_samples % 500 == 0)
            std::cout << number_of_samples << " samples to go" << std::endl;
      } while (number_of_samples);

      tee_log << "Max error found was: " << max_error << std::endl;

      unsigned current_bin = 0;
      unsigned lim = 1;
      unsigned old_lim = 0;

      while (errors.size())
      {
         old_lim = lim;
         lim *= 2;
         //std::cout << "Enter upper limit for bin " << current_bin << ": ";
         //std::cin >> lim;
         auto p = errors.upper_bound(error_data(0, 0, 0, lim));
         int bin_count = std::distance(errors.begin(), p);
         if (bin_count)
         {
            std::ofstream os((basename + "_errors_" + std::to_string(current_bin + 1) + ".csv").c_str());
            os << "a,b,z,error\n";
            bins[std::make_pair(old_lim, lim)] = bin_count;
            for (auto pos = errors.begin(); pos != p; ++pos)
            {
               os << pos->a << "," << pos->b << "," << pos->z << "," << pos->error << "\n";
            }
            errors.erase(errors.begin(), p);
         }
         ++current_bin;
      }

      tee_bins << "Results:\n\n";
      tee_bins << "#bin,Range,2^N,Count\n";
      int hash = 0;
      for (auto p = bins.begin(); p != bins.end(); ++p, ++hash)
      {
         tee_bins << hash << "," << p->first.first << "-" << p->first.second << "," << hash+1 << "," << p->second << std::endl;
      }
      if (evaluation_errors)
      {
         tee_bins << ",Failed,," << evaluation_errors << std::endl;
      }
      if (incalculable)
      {
         tee_bins << ",Incalculable,," << incalculable << std::endl;
      }
   }
   catch (const std::exception& e)
   {
      std::cout << "Terminating with unhandled exception: " << e.what() << std::endl;
   }

   return 0;
}

max_stars_repo_name: cpp-pm/boost

max_issues_repo_name: cpp-pm/boost

max_forks_repo_name: cpp-pm/boost


text:
#ifndef FILEFORMAT_PNG
#define FILEFORMAT_PNG

#include <png.h>
#include <zlib.h>

#include "image/format/base.hpp"
#include "image/image_base.hpp"

#include <boost/lexical_cast.hpp>

struct png_config {
  std::string comment;
  size_t center_x, center_y;
};

class png_format {
  public:
    typedef png_config opt_type;

    static void save(image_base* image, const std::string& path, opt_type& opts)
    {
      FILE *fp;
      png_structp write_struct = NULL;
      png_infop info_struct = NULL;
      png_bytep bytes_row = NULL;
      color* color_row = NULL;

      fp = fopen(path.c_str(), "wb");

      if (fp == NULL)
      {
        throw format_exception(strerror(errno));
      }

      write_struct = png_create_write_struct(PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);

      if (write_struct == NULL)
      {
        fclose(fp);
        throw format_exception("failed to acquire write_struct");
      }

      /* optipng tends to select these options for c10t generated images */
      png_set_compression_level(write_struct, 9);
      png_set_compression_mem_level(write_struct, 8);
      png_set_compression_strategy(write_struct, Z_DEFAULT_STRATEGY);
      png_set_filter(write_struct, PNG_FILTER_TYPE_BASE, PNG_FILTER_NONE);

      info_struct = png_create_info_struct(write_struct);

      if (info_struct == NULL)
      {
        fclose(fp);
        png_destroy_write_struct(&write_struct, NULL);
        throw format_exception("failed to acquire info_struct");
      }

      /* this is where we end up on errors */
      if (setjmp(png_jmpbuf(write_struct)))
      {
        fclose(fp);
        png_destroy_write_struct(&write_struct, &info_struct);
        delete [] bytes_row;
        delete [] color_row;
        throw format_exception("unknown libpng error");
      }

      png_init_io(write_struct, fp);

      png_set_IHDR(write_struct, info_struct, image->get_width(), image->get_height(),
            8, PNG_COLOR_TYPE_RGBA, PNG_INTERLACE_NONE,
            PNG_COMPRESSION_TYPE_BASE, PNG_FILTER_TYPE_BASE);

      if (!opts.comment.empty())
      {
         png_text title_text;
         title_text.compression = PNG_TEXT_COMPRESSION_NONE;
         title_text.key = (char *)"Title";
         title_text.text = (char *)opts.comment.c_str();
         png_set_text(write_struct, info_struct, &title_text, 1);
      }

      {
         png_text title_text;
         title_text.compression = PNG_TEXT_COMPRESSION_NONE;
         title_text.key = (char *)"center-x";
         std::string text = boost::lexical_cast<std::string, int>(opts.center_x);
         title_text.text = const_cast<char*>(text.c_str());
         png_set_text(write_struct, info_struct, &title_text, 1);
      }

      {
         png_text title_text;
         title_text.compression = PNG_TEXT_COMPRESSION_NONE;
         title_text.key = (char *)"center-y";
         std::string text = boost::lexical_cast<std::string, int>(opts.center_y);
         title_text.text = const_cast<char*>(text.c_str());
         png_set_text(write_struct, info_struct, &title_text, 1);
      }

      png_write_info(write_struct, info_struct);

      bytes_row = new png_byte[4 * image->get_width()];
      color_row = new color[image->get_width()];

      for (size_t y = 0; y < image->get_height(); y++)
      {
        image->get_line(y, color_row);

        for (unsigned int i = 0; i < image->get_width(); i++)
        {
          bytes_row[i*4 + 0] = png_byte(color_row[i].r * 255.0f);
          bytes_row[i*4 + 1] = png_byte(color_row[i].g * 255.0f);
          bytes_row[i*4 + 2] = png_byte(color_row[i].b * 255.0f);
          bytes_row[i*4 + 3] = png_byte(color_row[i].a * 255.0f);
        }

        png_write_row(write_struct, bytes_row);
      }

      png_write_end(write_struct, NULL);

      fclose(fp);
      png_free_data(write_struct, info_struct, PNG_FREE_ALL, -1);
      png_destroy_write_struct(&write_struct, &info_struct);
      delete [] bytes_row;
      delete [] color_row;
    }
};

#endif /* FILEFORMAT_PNG */