text:
#pragma once

#include "bytes_io_common.hpp"
#include <assert.h>
#include <boost/endian/conversion.hpp>
#include <cstring>
#include <stdint.h>
#include <type_traits>

namespace BytesIO {
template <ByteOrder byteOrder> class DynamicBytesWriter {
public:
  DynamicBytesWriter(size_t initialCapacity) : buffer_(new uint8_t[initialCapacity]), capacity_(initialCapacity), size_(0) {}

  ~DynamicBytesWriter() { delete[] buffer_; }

  void writeBytes(const uint8_t *bytes, size_t len) {
    ensureCapacity(len);
    memcpy(buffer_ + size_, bytes, len);
    size_ += len;
  }

  void writeChars(const char *chars, size_t len) {
    ensureCapacity(len);
    memcpy(buffer_ + size_, chars, len);
    size_ += len;
  }

  void setPosition(size_t position) { size_ = position; }

  template <class Enum, class = typename std::enable_if<std::is_enum<Enum>::value>::type> DynamicBytesWriter &operator<<(Enum u) { return *this << std::underlying_type_t<Enum>(u); }

  DynamicBytesWriter &operator<<(uint8_t u) { return write(u); }

  DynamicBytesWriter &operator<<(int8_t u) { return write(u); }

  DynamicBytesWriter &operator<<(char u) { return write(u); }

  DynamicBytesWriter &operator<<(uint16_t u) { return write(u); }

  DynamicBytesWriter &operator<<(int16_t u) { return write(u); }

  DynamicBytesWriter &operator<<(uint32_t u) { return write(u); }

  DynamicBytesWriter &operator<<(int32_t u) { return write(u); }

  DynamicBytesWriter &operator<<(uint64_t u) { return write(u); }

  DynamicBytesWriter &operator<<(int64_t u) { return write(u); }

  uint8_t *buffer() { return buffer_; }

  size_t size() { return size_; }

private:
  template <typename Integer> DynamicBytesWriter &write(Integer v) {
    ensureCapacity(sizeof(v));
    if (byteOrder == ByteOrder::BigEndian) {
      boost::endian::native_to_big_inplace(v);
    } else if (byteOrder == ByteOrder::LittleEndian) {
      boost::endian::native_to_little_inplace(v);
    }
    memcpy(buffer_ + size_, &v, sizeof(v));
    size_ += sizeof(v);
    return *this;
  }

  void unsafeWriteU8(uint8_t u) { buffer_[size_++] = u; }

  void ensureCapacity(size_t bytesToWrite) { ensureTotalCapacity(size_ + bytesToWrite); }

  void ensureTotalCapacity(size_t totalCapacity) {
    if (totalCapacity > capacity_) {
      size_t targetCapacity = capacity_;
      while (totalCapacity > targetCapacity) {
        // ~1.5x
        targetCapacity = (targetCapacity * 3 + 2) / 2;
      }
      enlarge(targetCapacity);
    }
  }

  void enlarge(size_t newCapacity) {
    uint8_t *newBuffer = new uint8_t[newCapacity];
    std::memmove(newBuffer, buffer_, capacity_);
    buffer_ = newBuffer;
    capacity_ = newCapacity;
  }

  uint8_t *buffer_;
  size_t size_;
  size_t capacity_;
};
} // namespace BytesIO


text:
#ifndef OPENP2P_CONDITION_HPP
#define OPENP2P_CONDITION_HPP

#include <boost/thread.hpp>

#include <OpenP2P/Lock.hpp>
#include <OpenP2P/Timeout.hpp>

namespace OpenP2P {

	class Condition {
		public:
			inline void notifyOne() {
				internalCondition_.notify_one();
			}

			inline void notifyAll() {
				internalCondition_.notify_all();
			}

			inline bool wait(Lock& lock, Timeout timeout = Timeout::Infinite()) {
				if (timeout.isInfinite()) {
					internalCondition_.wait(lock);
					return true;
				} else {
					return internalCondition_.timed_wait(lock, boost::posix_time::milliseconds(timeout.seconds() * 1000));
				}
			}

		private:
			boost::condition_variable internalCondition_;

	};

}

#endif

ext: hpp


text:
/*
 * cte.cpp
 *
 *  Created on: 14 Sep 2018
 *      Author: radu
 *
 * Copyright (c) 2015, International Business Machines Corporation
 * and University of California Irvine. All rights reserved.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
 * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
 * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
 * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */

#include "cte.h"

#include <boost/graph/adjacency_list.hpp>
#include <boost/graph/kruskal_min_spanning_tree.hpp>

using namespace boost;
typedef adjacency_list<vecS, vecS, undirectedS, no_property,
		property<edge_weight_t, int> > BoostGraph;
typedef graph_traits<BoostGraph>::edge_descriptor BoostEdge;
typedef graph_traits<BoostGraph>::vertex_descriptor BoostVertex;
typedef std::pair<int, int> E;

namespace merlin {

// Initialize the clique tree
void cte::init() {

	// Prologue
	if (m_verbose > 0) {
		std::cout << "[CTE] + inference task   : " << m_task << std::endl;
	}

	if (m_query.empty() == false) {
		if (m_verbose > 0) {
			std::cout << "[CTE] + query vars       : ";
			std::copy(m_query.begin(), m_query.end(), std::ostream_iterator<vindex>(std::cout, " "));
			std::cout << std::endl;
		}

		// Add a dummy factor over the query variables
		if (m_query.size() <= MERLIN_MAXSIZE_JOINT_MARGINAL) {
			variable_set scope;
			for (size_t i = 0; i < m_query.size(); ++i) {
				scope |= m_gmo.var(m_query[i]);
			}

			factor q(scope, 1.0);
			m_gmo.add_factor(q);
		}
	}

	if (m_order.size() == 0) { // if we need to construct an elimination ordering
		m_order = m_gmo.order(m_order_method);
	}

	if (m_verbose > 0) {
		std::cout << "[CTE] + ordering method  : " << m_order_method << std::endl;
		std::cout << "[CTE] + elimination      : ";
		std::copy(m_order.begin(), m_order.end(),
				std::ostream_iterator<size_t>(std::cout, " "));
		std::cout << std::endl;

		// Calculate the induced width of the elimination ordering
		size_t wstar = m_gmo.induced_width(m_order);
		std::cout << "[CTE] + induced width    : " << wstar << std::endl;
		std::cout << "[CTE] + exact inference  : Yes" << std::endl;
		std::cout << "[CTE] + ordering time    : " << (timeSystem() - m_start_time) << " seconds" << std::endl;
		std::cout << "[CTE] Building clique tree ... " << std::endl;
	}

	// Get the factors scopes
	size_t n = m_gmo.num_nodes();
	vector<variable_set> fin;
	for (vector<factor>::const_iterator i = m_gmo.get_factors().begin();
			i != m_gmo.get_factors().end(); ++i) {
		fin.push_back((*i).vars());
	}

	// Create the graph
	graph g(n);
	g.init(fin);

	// Triangulate the graph
	g.triangulate(m_order);

	// Find the maximal cliques
	std::vector<std::set<size_t> > clusters = g.maximal_cliques(m_order);

	// Build the clique tree
	build_clique_tree(clusters);

	if (m_verbose > 0) {
		std::cout << "[CTE] Created clique tree with " << m_clusters.size()
				<< " clique factors" << std::endl;
	}

	size_t max_clique_size = 0, max_sepset_size = 0;
	for (size_t i = 0; i < m_edges.size(); ++i) {
		detail::edge* e = m_edges[i];
		max_clique_size = std::max(max_clique_size, e->first->clique.size());
		max_clique_size = std::max(max_clique_size, e->second->clique.size());
		max_sepset_size = std::max(max_sepset_size, e->sepset.size());
	}

	// Initialize beliefs (marginals)
	m_logz = 0;
	m_beliefs.clear();
	m_beliefs.resize(m_gmo.nvar(), factor(1.0));

	// Output the bucket tree statistics
	if (m_verbose > 0) {
		double elapsed = (timeSystem() - m_start_time);
		std::cout << "[CTE] Number of cliques  : " << m_clusters.size() << std::endl;
		std::cout << "[CTE] Number of edges    : " << m_edges.size() << std::endl;
		std::cout << "[CTE] Max clique size    : " << max_clique_size << std::endl;
		std::cout << "[CTE] Max separator size : " << max_sepset_size << std::endl;
		std::cout << "[CTE] Finished initialization in " << elapsed << " seconds" << std::endl;
	}
}

// Build the join tree and initialize the cluster factors
void cte::build_clique_tree(std::vector<std::set<size_t> >& cliques) {

	size_t n = cliques.size();
	std::vector<std::vector<int> > adjacencies;
	adjacencies.resize(n);
	for (std::vector<std::vector<int> >::iterator it = adjacencies.begin();
			it != adjacencies.end(); ++it) {
		it->resize(n);
	}

	// create nodes
	m_clusters.resize(n);
	for (size_t i = 0; i < n; ++i) {
		m_clusters[i].id = i;
		for (std::set<size_t>::iterator j = cliques[i].begin();
				j != cliques[i].end(); ++j) {
			m_clusters[i].clique |= var(*j);
		}
	}

	// fill in the weights
	for (size_t i = 0; i < n; ++i) {
		for (size_t j = 0; j < n; ++j) {
			if (i != j) {
				detail::node& ni = m_clusters[i];
				detail::node& nj = m_clusters[j];

				std::set<int> sep;
				std::set_intersection(ni.clique.begin(), ni.clique.end(),
						nj.clique.begin(), nj.clique.end(),
						std::inserter(sep, sep.begin()));

				if (sep.empty()) {
					adjacencies[i][j] = (100000);
				} else {
					adjacencies[i][j] = -(sep.size());
				}
			}
		}
	}

	// create a boost weighted graph and run Kruskal on it
	const int num_nodes = n;
	std::vector<E> temp1;
	std::vector<int> temp2;
	for (int i = 0; i < num_nodes - 1; ++i) {
		for (int j = i + 1; j < num_nodes; ++j) {
			if (i != j) {
				temp1.push_back( E(i, j) );
				temp2.push_back( adjacencies[i][j] );
			}
		}
	}

	E* edge_array = new E[temp1.size()];
	int* weights = new int[temp1.size()];
	for (size_t i = 0; i < temp1.size(); ++i) {
		edge_array[i] = temp1[i];
		weights[i] = temp2[i];
	}

	std::size_t num_edges = temp1.size();

	BoostGraph g(edge_array, edge_array + num_edges, weights, num_nodes);
	std::vector<BoostEdge> spanning_tree;

	kruskal_minimum_spanning_tree(g, std::back_inserter(spanning_tree));

	if (m_debug) {
		std::cout << "Undirected junction tree (MST):" << std::endl;
		property_map<BoostGraph, edge_weight_t>::type weight = get(edge_weight, g);
		for (std::vector<BoostEdge>::iterator ei = spanning_tree.begin();
				ei != spanning_tree.end(); ++ei) {
			std::cout << source(*ei, g) << " <--> " << target(*ei, g)
					<< " with weight of " << weight[*ei] << "\n";
		}
	}

	// select the root and redirect edges outwards
	m_root = &(m_clusters.back());
	size_t root = m_root->id;
	std::set<size_t> visited;
	std::stack<size_t> dfs;
	dfs.push(root);
	while ( !dfs.empty() ) {

		size_t c = dfs.top(); // current node
		dfs.pop();

		// get all unvisited children and direct them towards n
		for (std::vector<BoostEdge>::iterator ei = spanning_tree.begin();
				ei != spanning_tree.end(); ++ei) {
			size_t src = (size_t) source(*ei, g);
			size_t trg = (size_t) target(*ei, g);
			detail::node* from = NULL, *to = NULL;

			if (src == c &&
				std::find(visited.begin(), visited.end(), trg) == visited.end()) {

				from = &(m_clusters[trg]); // from
				to = &(m_clusters[src]); // to
				dfs.push(trg);
			}

			if (trg == c &&
				std::find(visited.begin(), visited.end(), src) == visited.end()) {

				from = &(m_clusters[src]); // from
				to = &(m_clusters[trg]); // to
				dfs.push(src);
			}

			// create the directed edge
			if (from != NULL && to != NULL) {
				detail::edge* e = new detail::edge(from, to);
				m_edges.push_back(e);
				from->edges.push_back(e);
				to->edges.push_back(e);
				to->children.push_back(from);
				from->parent = to;
			}
		}

		visited.insert(c);
	}

	// find the message order (schedule) by a bfs of the junction tree
	std::queue<detail::node*> bfs;
	bfs.push(m_root);
	while (!bfs.empty()) {
		detail::node* c = bfs.front();
		bfs.pop();
		vector<detail::edge*>& elist = c->edges;
		for (vector<detail::edge*>::iterator it = elist.begin();
				it != elist.end(); ++it) {
			detail::edge* e = (*it);
			if (e->second->id == c->id) {
				m_messages.push_back(e);
			}
		}

		for (vector<detail::node*>::iterator it = c->children.begin();
				it != c->children.end(); ++it) {
			bfs.push(*it);
		}
	}

	// safety checks
	assert(m_messages.size() == m_edges.size());

	// reverse the order to start from the leaves
	std::reverse(m_messages.begin(), m_messages.end());

	// clean up
	delete[] weights;
	delete[] edge_array;

	if (m_debug) {
		std::ofstream fout("kruskal.dot");
		fout << "digraph JT {\n" << " rankdir=LR\n" << " size=\"3,3\"\n"
				<< " ratio=\"filled\"\n" << " edge[style=\"bold\"]\n"
				<< " node[shape=\"circle\"]\n";
		// nodes
		for (size_t i = 0; i < m_clusters.size(); ++i) {
			std::stringstream label;
			label << m_clusters[i].id << ": " << m_clusters[i].clique;
			fout << "node" << m_clusters[i].id
				<< "[ label = \"" << label.str() << "\"];\n";
		}

		// edges
		for (size_t i = 0; i < m_edges.size(); ++i) {
			detail::edge* e = m_edges[i];
			fout << "node" << e->first->id << " -> " << "node" << e->second->id << ";\n";
		}

		fout << "}\n";
		fout.close();
	}

	// Map variables to clusters
	m_var2clique.resize(m_gmo.nvar(), -1);
	for (size_t i = 0; i < m_clusters.size(); ++i) {
		detail::node& cl = m_clusters[i];
		variable_set& vars = cl.clique;
		for (variable_set::const_iterator j = vars.begin();
				j != vars.end(); ++j) {
			size_t v = (*j).label();
			if (m_var2clique[v] < 0) {
				m_var2clique[v] = cl.id;
			}
		}
	}

	// Allocate original functions to clusters and compute clique factors
	size_t idx = 0;
	for (vector<factor>::const_iterator fi = m_gmo.get_factors().begin();
			fi != m_gmo.get_factors().end(); ++fi, ++idx) {

		const factor& f = (*fi);
		for (size_t ci = 0; ci != m_clusters.size(); ++ci) {
			detail::node& cl = m_clusters[ci];
			if ((f.vars() & cl.clique) == f.vars()) {
				cl.originals.push_back(idx);
				cl.theta *= f;
				break;
			}
		}
	}

	if (m_debug) {
		std::cout << "Initial cluster factors:" << std::endl;
		for (size_t ci = 0; ci != m_clusters.size(); ++ci) {
			detail::node& cl = m_clusters[ci];
			std::cout << " " << ci << " -> " << cl.theta << std::endl;
		}
	}
}

// Update the collection of factors, recompute the clique factors and reset the
// messages. We assume there is a one-to-one mapping between the old factors
// and the new factors. This method is to be used only for EM param learning.
void cte::reinit(const std::vector<factor>& factors) {

	// Set the new factors
	for (size_t i = 0; i < factors.size(); ++i) {
		m_gmo.set_factor(i, factors[i]);
	}

	// Recompute the clique factors
	for (size_t ci = 0; ci != m_clusters.size(); ++ci) {
		detail::node& cl = m_clusters[ci];
		cl.theta = factor(1.0);
		for (size_t j = 0; j < cl.originals.size(); ++j) {
			const factor& f = m_gmo.get_factor(cl.originals[j]);
			cl.theta *= f;
		}
	}

	// Reset the edge messages
	for (vector<detail::edge*>::iterator it = m_edges.begin();
			it != m_edges.end(); ++it) {
		(*it)->reset();
	}

	// Reset the log partition function and beliefs
	m_logz = 0;
	m_beliefs.clear();
	m_beliefs.resize(m_gmo.nvar(), factor(1.0));

	if (m_debug) {
		std::cout << "Initial cluster factors:" << std::endl;
		for (size_t ci = 0; ci != m_clusters.size(); ++ci) {
			detail::node& cl = m_clusters[ci];
			std::cout << " " << ci << " -> " << cl.theta << std::endl;
		}
	}
}

// Forward message propagation (from leaves to root)
void cte::forward() {
	double timestamp = timeSystem();
	for (vector<detail::edge*>::iterator i = m_messages.begin();
			i != m_messages.end(); ++i) {
		detail::edge* e = (*i);
		e->messageFwd();
	}

	if (m_verbose > 0) {
		std::cout << "[CTE] Finished forward pass in " << (timeSystem() - timestamp)
			<< " seconds" << std::endl;
	}
}

// Forward message propagation with evidence (from leaves to root)
void cte::forward(std::vector<int>& evidence) {
	for (vector<detail::edge*>::iterator i = m_messages.begin();
			i != m_messages.end(); ++i) {
		detail::edge* e = (*i);
		e->messageFwd(evidence);

		if (m_debug) {
			std::cout << " -> forward msg from "
				<< e->first->id << " to " << e->second->id
				<< ": " << e->fwd << std::endl;
		}
	}
}

// Forward message on an edge
void detail::edge::messageFwd() {

	// collect the original functions
	factor F = first->theta;

	// collect the incoming messages (should already be projected)
	for (size_t i = 0; i < first->edges.size(); ++i) {
		detail::edge* e = first->edges[i];
		if (e == this) {
			continue; // skip current edge
		}
		if (e->first->id == first->id) {
			F *= e->bwd; //e->getMessage2();
		} else if (e->second->id == first->id) {

			F *= e->fwd; //e->getMessage1();
		}
	}

	// marginalize the eliminator
	variable_set elim = (first->clique - sepset);
	if (elim.size() > 0) {
		F = F.sum(elim);
	}

	// update the forward message
	this->fwd = F;
}

// Forward message computed with evidence on an edge
void detail::edge::messageFwd(std::vector<int>& evidence) {

	// collect the original functions
	factor F = first->theta;
	F = F.condition(evidence);

	// collect the incoming messages (should already be projected)
	for (size_t i = 0; i < first->edges.size(); ++i) {
		detail::edge* e = first->edges[i];
		if (e == this) {
			continue; // skip current edge
		}
		if (e->first->id == first->id) {
			F *= e->bwd.condition(evidence);
		} else if (e->second->id == first->id) {

			F *= e->fwd.condition(evidence);
		}
	}

	// Find the evidence variables in the scope of the cluster
	variable_set evid;
	for (variable_set::const_iterator i = first->clique.begin();
			i != first->clique.end(); ++i) {
		int v = i->label();
		if (evidence[v] >= 0) {
			evid |= *i;
		}
	}

	// Marginalize the non-evidence eliminator
	//F = F.condition(evidence);
	variable_set elim = (first->clique - sepset);
	elim -= evid;
	if (elim.size() > 0) {
		F = F.sum(elim);
	}

	// update the forward message
	this->fwd = F;

}

// Backward message propagation (from root to leaves)
void cte::backward() {
	double timestamp = timeSystem();
	for (vector<detail::edge*>::reverse_iterator ri = m_messages.rbegin();
			ri != m_messages.rend(); ++ri) {
		detail::edge* e = (*ri);
		e->messageBwd();
	}

	if (m_verbose > 0) {
		std::cout << "[CTE] Finished backward pass in " << (timeSystem() - timestamp)
			<< " seconds" << std::endl;
	}
}

// Backward message propagation with evidence (from root to leaves)
void cte::backward(std::vector<int>& evidence) {
	for (vector<detail::edge*>::reverse_iterator ri = m_messages.rbegin();
			ri != m_messages.rend(); ++ri) {
		detail::edge* e = (*ri);
		e->messageBwd(evidence);

		if (m_debug) {
			std::cout << " <- backward msg from "
				<< e->second->id << " to " << e->first->id
				<< ": " << e->bwd << std::endl;
		}
	}
}

// Backward message on an edge
void detail::edge::messageBwd() {

	// collect the original functions
	factor F = second->theta;

	// collect the incoming messages
	for (size_t i = 0; i < second->edges.size(); ++i) {
		detail::edge* e = second->edges[i];
		if (e == this) {
			continue; // skip current edge
		}
		if (e->first->id == second->id) {
			F *= e->bwd; //e->getMessage2();
		} else if (e->second->id == second->id) {
			F *= e->fwd; //e->getMessage1();
		}
	}

	// marginalize the eliminator
	variable_set elim = (second->clique - sepset);
	if (elim.size() > 0) {
		F = F.sum(elim);
	}

	this->bwd = F;
}

// Backward message computed with evidence on an edge
void detail::edge::messageBwd(std::vector<int>& evidence) {

	// collect the original functions
	factor F = second->theta;
	F = F.condition(evidence);

	// collect the incoming messages
	for (size_t i = 0; i < second->edges.size(); ++i) {
		detail::edge* e = second->edges[i];
		if (e == this) {
			continue; // skip current edge
		}
		if (e->first->id == second->id) {
			F *= e->bwd.condition(evidence); //e->getMessage2();
		} else if (e->second->id == second->id) {
			F *= e->fwd.condition(evidence); //e->getMessage1();
		}
	}

	// find the evidence variables in the scope of the cluster
	variable_set evid;
	for (variable_set::const_iterator i = second->clique.begin();
			i != second->clique.end(); ++i) {
		int v = i->label();
		if (evidence[v] >= 0) {
			evid |= *i;
		}
	}

	// marginalize the eliminator
	variable_set elim = (second->clique - sepset);
	elim -= evid;
	if (elim.size() > 0) {
		F = F.sum(elim);
	}

	this->bwd = F;

}

// Clique tree calibration
void cte::calibrate() {

	forward();
	backward();
}

// Update the beliefs (marginals or max-marginals)
void cte::update() {

	// Compute the clique beliefs
	for (size_t i = 0; i < m_clusters.size(); ++i) {
		detail::node& cl = m_clusters[i];
		cl.belief = cl.theta;
		for (vector<detail::edge*>::iterator ei = cl.edges.begin();
				ei != cl.edges.end(); ++ei) {
			detail::edge* e = (*ei);
			if (e->second->id == cl.id) {
				cl.belief *= e->fwd;
			}
			if (e->first->id == cl.id) {
				cl.belief *= e->bwd;
			}
		}
	}

	// Compute the log partition function
	m_logz = std::log(m_root->belief.sum());

	// Compute the marginals
	for (size_t v = 0; v < m_var2clique.size(); ++v) {
		size_t ci = m_var2clique[v];
		detail::node& cl = m_clusters[ci];
		variable VX = m_gmo.var(v);

		m_beliefs[v] = marg(cl.belief, VX);
		m_beliefs[v].normalize(); // normalize
	}
}

// Compute the joint marginal of a set of variables
void cte::joint_marginal(const variable_set& scope) {

	// Find the shallowest cliques that contain the scope
	std::vector<int> nodes;
	variable_set temp = scope;
	while (temp.size() > 0) {
		int max_id = -1, max_score = -1;
		std::queue<detail::node*> bfs;
		bfs.push(m_root);
		while (!bfs.empty()) {
			detail::node* c = bfs.front();
			bfs.pop();

			int score = (c->clique & temp).size();
			if (score > 0 && score > max_score) {
				max_score = score;
				max_id = c->id;
			}

			for (vector<detail::node*>::iterator it = c->children.begin();
					it != c->children.end(); ++it) {
				bfs.push(*it);
			}
		}

		assert(max_id >= 0);
		nodes.push_back(max_id);
		temp = temp - (temp & m_clusters[max_id].clique);
	}

	if (m_debug) {
		std::cout << "[DEBUG] Found the shallowest nodes: ";
		std::copy(nodes.begin(), nodes.end(),
				std::ostream_iterator<int>(std::cout, " "));
		std::cout << std::endl;
	}

	// Check if the scope is included in one clique
	if (nodes.size() == 1) { // one clique
		m_marginal = marg(m_clusters[nodes.at(0)].belief, scope);
		m_marginal.normalize();
	} else { // multiple cliques

		// Collect all relevant factors
		vector<factor> factors;
		factors.push_back(m_root->belief);
		for (size_t i = 0; i < nodes.size(); ++i) {
			detail::node* c = &(m_clusters.at(nodes[i]));
			while (c != NULL) {
				for (vector<detail::edge*>::iterator ei = c->edges.begin();
						ei != c->edges.end(); ++ei) {
					detail::edge* e = (*ei);
					if (e->first->id == c->id) {
						factor f = c->belief;
						f /= e->fwd;
						factors.push_back(f);
						break;
					}
				}

				c = c->parent;
				if (c == m_root) {
					break;
				}
			}
		}

		variable_set all_vars, elim_vars;
		for (vector<factor>::iterator fi = factors.begin();
				fi != factors.end(); ++fi) {
			all_vars |= fi->vars();
		}
		elim_vars = (all_vars - scope);

		// Run variables elimination over these factors
		variable_order_t order;
		for (size_t i = 0; i < m_order.size(); ++i) {
			size_t v = m_order[i];
			if (elim_vars.contains(m_gmo.var(v))) {
				order.push_back(v);
			}
		}

		if (m_debug) {
			std::cout << "[DEBUG] All vars: " << all_vars << std::endl;
			std::cout << "[DEBUG] Scope: " << scope << std::endl;
			std::cout << "[DEBUG] Elim: " << elim_vars << std::endl;
			std::cout << "[DEBUG] Factors: " << factors.size() << std::endl;
			std::cout << "[DEBUG] Elim order: ";
			std::copy(order.begin(), order.end(),
					std::ostream_iterator<size_t>(std::cout, " "));
			std::cout << std::endl;
		}

		// Run variable elimination
		for (size_t i = 0; i < order.size(); ++i) {
			size_t v = order[i];
			variable VX = m_gmo.var(v);

			// Collect all factors mentioning VX
			factor f(1.0);
			for (vector<factor>::iterator it = factors.begin();
					it != factors.end();) {
				if (it->variables().contains(VX)) {
					f *= (*it);
					it = factors.erase(it);
				} else {
					++it;
				}
			}

			// Eliminate the variable
			f = elim(f, VX);
			factors.push_back(f);
		}

		// Compute the joint marginal
		m_marginal = factor(1.0);
		for (vector<factor>::iterator it = factors.begin();
				it != factors.end(); ++it) {
			m_marginal *= (*it);
			m_marginal.normalize();
		}
	}
}

// The scope must be included in one cluster (Bayes nets only)
void cte::joint_marginal(const variable_set& scope, std::vector<int>& evidence) {

	// Get the cluster that contains the scope
	bool found = false;
	size_t j = 0;
	for (size_t i = 0; i < m_clusters.size(); ++i) {
		detail::node& cl = m_clusters[i];
		variable_set temp = (cl.clique & scope);
		if (temp == scope) {
			found = true;
			j = i;
			break;
		}
	}

	// Safety check
	assert(found);

	// Compute the joint marginal (belief) subject to evidence
	m_marginal = factor(scope, 0.0);
	detail::node& n = m_clusters[j];
	n.belief = n.theta.condition(evidence);
	for (vector<detail::edge*>::iterator ei = n.edges.begin();
			ei != n.edges.end(); ++ei) {
		detail::edge* e = (*ei);
		if (e->second->id == n.id) {
			n.belief *= e->fwd;
		}
		if (e->first->id == n.id) {
			n.belief *= e->bwd;
		}
	}

	// Update the belief entries subject to evidence
	variable_set v_rem;
	for (variable_set::const_iterator vi = n.clique.begin();
			vi != n.clique.end(); ++vi) {
		if (!scope.contains(*vi) && evidence[vi->label()] < 0) {
			v_rem |= *vi;
		}
	}
	n.belief = n.belief.sum(v_rem);

	if (m_debug) {
		std::cout << "[DEBUG] Joint marginal scope: " << m_marginal.vars() << std::endl;
		std::cout << "[DEBUG] Actual belief scope:  " << n.belief << std::endl;
	}

	// Fill in all configurations of the marginal (including the evidence)
	index_config cv1(m_marginal.vars(), true);
	config_index cv2(n.belief.vars(), true);
	for (size_t i = 0; i < m_marginal.num_states(); ++i) {
		std::map<size_t, size_t> config = cv1.convert(i);
		if (is_compatible(config, evidence)) {
			size_t j = cv2.convert(config);
			double v = n.belief.get(j);
			m_marginal.set(i, v);
		} else {
			m_marginal.set(i, 0.0);
		}
	}

	// Normalize by P(evidence)
	m_marginal /= std::exp(m_logz);

	if (m_debug) {
		std::cout << "[DEBUG] Joint marginal: " << m_marginal << std::endl;
	}
}

// Check if a variable configuration is compatible with the evidence vector
bool cte::is_compatible(std::map<size_t, size_t>& config,
		std::vector<int>& evidence) {
	bool ok = true;
	for (std::map<size_t, size_t>::iterator mi = config.begin();
			mi != config.end(); ++mi) {
		size_t var = mi->first;
		size_t val = mi->second;
		if (evidence[var] < 0) {
			continue; // skip non-evidence variable
		} else {
			if (val != evidence[var]) {
				ok = false;
				break;
			}
		}
	}

	return ok;
}

// Run the clique-tree algorithm with evidence (for EM learning)
// Returns *true* if P(evidence) > 0, and *false* otherwise
bool cte::propagate_evidence(std::vector<int>& evidence) {

	// Save the current evidence
	if (m_debug) {
		std::cout << "[CTE] Propagate evidence: ";
		std::copy(evidence.begin(), evidence.end(), std::ostream_iterator<int>(std::cout, " "));
		std::cout << std::endl;
	}

	m_evidence = evidence;
	forward(evidence);
	backward(evidence);

	// Compute the posterior marginals for all variables (evidence, non-evidence)
	for (size_t i = 0; i < m_clusters.size(); ++i) {
		detail::node& cl = m_clusters[i];
		cl.belief = cl.theta.condition(evidence);
		for (vector<detail::edge*>::iterator ei = cl.edges.begin();
				ei != cl.edges.end(); ++ei) {
			detail::edge* e = (*ei);
			if (e->second->id == cl.id) {
				cl.belief *= e->fwd;
			}
			if (e->first->id == cl.id) {
				cl.belief *= e->bwd;
			}
		}
	}

	// Compute the root clique belief (to get log partition function)
	m_root->belief = m_root->theta.condition(evidence);
	for (vector<detail::edge*>::iterator ei = m_root->edges.begin();
		ei != m_root->edges.end(); ++ei) {
		detail::edge* e = (*ei);
		if (e->second->id == m_root->id) {
			m_root->belief *= e->fwd;
		}
		if (e->first->id == m_root->id) {
			m_root->belief *= e->bwd;
		}
	}

	// Compute the probability of evidence
	bool result = true;
	factor::value pe = m_root->belief.sum();
	if (pe == 0.0) {
		result = false;
	}

	// Update the log partition function
	m_logz = std::log(m_root->belief.sum());

	// Compute the marginals
	for (size_t v = 0; v < m_var2clique.size(); ++v) {
		size_t ci = m_var2clique[v];
		detail::node& cl = m_clusters[ci];
		variable VX = m_gmo.var(v);

		if (evidence[v] >= 0) {
			factor tmp(variable_set(VX), 0.0);
			tmp.set(evidence[v], 1.0);
			m_beliefs[v] = tmp;
		} else {
			m_beliefs[v] = marg(cl.belief, VX);
			m_beliefs[v].normalize(); // normalize
		}
	}

	if (m_debug) {
		std::cout << "[CTE] Finished propagating evidence with logZ = "
				<< m_logz << " (" << std::exp(m_logz) << ")" << std::endl;
		std::cout << "[CTE] Posterior marginals:" << std::endl;
		for (size_t i = 0; i < m_beliefs.size(); ++i) {
			std::cout << " " << m_beliefs[i] << std::endl;
		}
	}

	return result;
}

/// Run the clique-tree elimination algorithm.
void cte::run() {

	// Start the timer and store it
	m_start_time = timeSystem();

	init();
	calibrate();
	update();
	if (!m_query.empty()) {
		variable_set scope;
		for (size_t i = 0; i < m_query.size(); ++i) {
			scope |= m_gmo.var(m_query[i]);
		}

		joint_marginal(scope);
	}

	// Output solution (UAI output format)
	if (m_verbose > 0) {
		std::cout << "[CTE] Finished in " << (timeSystem() - m_start_time)
			<< " seconds" << std::endl;
	}

	// Compute the joint marginal given in query

	switch (m_task) {
	case Task::PR:	// partition function
		{
			std::cout << "PR" << std::endl;
			std::cout << std::fixed << std::setprecision(MERLIN_PRECISION)
				<< m_logz << " (" << std::scientific
				<< std::setprecision(MERLIN_PRECISION)
				<< std::exp(m_logz) << ")" << std::endl;

			if (isinf(m_logz)) {
				std::cout << "STATUS" << std::endl;
				std::cout << "false: Inconsistent evidence or underflow" << std::endl;
			} else {
				std::cout << "STATUS" << std::endl;
				std::cout << "true: Consistent evidence" << std::endl;
			}

			break;
		}
	case Task::MAR:	// marginals
		{
			std::cout << "PR" << std::endl;
			std::cout << std::fixed << std::setprecision(MERLIN_PRECISION)
				<< m_logz << " (" << std::scientific
				<< std::setprecision(MERLIN_PRECISION)
				<< std::exp(m_logz) << ")" << std::endl;

			if (isinf(m_logz)) {
				std::cout << "STATUS" << std::endl;
				std::cout << "false: Inconsistent evidence or underflow" << std::endl;
			} else {
				std::cout << "STATUS" << std::endl;
				std::cout << "true: Consistent evidence" << std::endl;
			}

			std::cout << "MAR" << std::endl;
			std::cout << m_gmo.nvar();
			for (vindex v = 0; v < m_gmo.nvar(); ++v) {
				variable VX = m_gmo.var(v);
				std::cout << " " << VX.states();
				for (size_t k = 0; k < VX.states(); ++k) {
					std::cout << " " << std::fixed
						<< std::setprecision(MERLIN_PRECISION)
						<< belief(VX)[k];
				}
			}
			std::cout << std::endl;
			if (!m_query.empty()) {
				variable_set scope = m_marginal.variables();
				std::cout << "JOINT_MAR : " << scope << std::endl;
				std::vector<size_t> dims(scope.size());
				size_t j = 0;
				for (variable_set::const_iterator si = scope.begin(); si != scope.end(); ++si) {
					dims[j++] = (*si).states();
				}

				for (size_t index = 0; index < m_marginal.numel(); ++index) {
					std::vector<size_t> I(dims.size(), 0); // configuration of source variables
					size_t i = index;
					for (size_t v = 0; v < dims.size(); ++v) {
						I[v] = i % dims[v];
						i -= I[v];
						i /= dims[v];
					}

					for (size_t j = 0; j < I.size(); ++j) {
						std::cout << I[j] << " ";
					}
					std::cout << ": " << std::fixed
							<< std::setprecision(MERLIN_PRECISION)
							<< m_marginal[index]
							<< std::endl;
				}
			}

			break;
		}
	default:
		break;
	}
}

/// Write the solution to the output stream
void cte::write_solution(std::ostream& out, const std::map<size_t, size_t>& evidence,
		const std::map<size_t, size_t>& old2new, const graphical_model& orig,
		const std::set<size_t>& dummies, int output_format) {

	if (output_format == MERLIN_OUTPUT_JSON) {
		out << "{";
		out << " \"algorithm\" : \"cte\", ";
		switch (m_task) {
		case Task::PR:
			{
				double val = m_logz + std::log(orig.get_global_const());
				double prob = std::exp(val);

				out << " \"task\" : \"PR\", ";
				out << " \"value\" : " << std::fixed
					<< std::setprecision(MERLIN_PRECISION)
					<< (m_logz + std::log(orig.get_global_const())) << ", ";
				if (prob == 0.0) { // probability of evidence is 0
					out << " \"status\" : \"false\", ";
					out << " \"message\" : \"Inconsistent evidence or underflow\" ";
				} else {
					out << " \"status\" : \"true\", ";
					out << " \"message\" : \"Consistent evidence\" ";
				}

				break;
			}
		case Task::MAR:
			{
				double val = m_logz + std::log(orig.get_global_const());
				double prob = std::exp(val);

				out << " \"task\" : \"MAR\", ";
				out << " \"value\" : " << std::fixed
					<< std::setprecision(MERLIN_PRECISION)
					<< (m_logz + std::log(orig.get_global_const())) << ", ";

				if (prob == 0.0) { // probability of evidence is 0
					out << " \"status\" : \"false\", ";
					out << " \"message\" : \"Inconsistent evidence or underflow\", ";
					out << " \"marginals\" : [] ";
				} else {
					out << " \"status\" : \"true\", ";
					out << " \"message\" : \"Consistent evidence\", ";
					out << " \"marginals\" : [ ";

					for (vindex i = 0; i < orig.nvar(); ++i) {
						if (dummies.find(i) != dummies.end()) {
							continue; // skip dummy variables from output
						}

						variable v = orig.var(i);
						out << "{";
						out << " \"variable\" : " << v.label() << ", ";
						out << " \"states\" : " << v.states() << ", ";
						out << " \"probabilities\" : [";
						try { // evidence variable
							size_t val = evidence.at(i);
							for (size_t k = 0; k < v.states(); ++k) {
								out << std::fixed
									<< std::setprecision(MERLIN_PRECISION)
									<< (k == val ? 1.0 : 0.0);
								if (k != v.states() - 1) {
									out << ", ";
								}
							}
							out << "] ";
						} catch(std::out_of_range& e) { // non-evidence variable
							vindex vx = old2new.at(i);
							variable VX = var(vx);
							for (size_t k = 0; k < VX.states(); ++k) {
								out << std::fixed
									<< std::setprecision(MERLIN_PRECISION)
									<< belief(VX)[k];
								if (k != VX.states() - 1) {
									out << ", ";
								}
							}
							out << "] ";
						}
						out << "}";
						if (i != orig.nvar() - 1) {
							out << ", ";
						}
					} // end for
					out << "] ";

					if (!m_query.empty()) {
						out << ", ";
						out << "\"joint_marginal\" : {";
						variable_set scope = m_marginal.variables();
						out << "\"scope\" : [";
						for (size_t i = 0; i < m_query.size(); ++i) {
							out << m_query[i];
							if (i < m_query.size() - 1) {
								out << ",";
							}
						}
						out << "], " << std::endl;
						out << "\"probabilities\" : [";

						std::vector<size_t> dims(scope.size());
						size_t j = 0;
						for (variable_set::const_iterator si = scope.begin();
								si != scope.end(); ++si) {
							dims[j++] = (*si).states();
						}

						size_t N = m_marginal.numel();
						for (size_t index = 0; index < N; ++index) {
							std::vector<size_t> I(dims.size(), 0); // configuration of source variables
							size_t i = index;
							for (size_t v = 0; v < dims.size(); ++v) {
								I[v] = i % dims[v];
								i -= I[v];
								i /= dims[v];
							}

							out << "{\"config\" : [";
							for (size_t j = 0; j < I.size(); ++j) {
								out << I[j];
								if (j < I.size() - 1) {
									out << ",";
								}
							}
							out << "], ";
							out << "\"value\" : " << std::fixed
									<< std::setprecision(MERLIN_PRECISION)
									<< m_marginal[index]
									<< "}";
							if (index < N - 1) {
								out << ", ";
							}
						}
						out << "]}";
					}
				}

				break;
			}
		default:
			break;
		}

		out << "}";
	} else if (output_format == MERLIN_OUTPUT_UAI) {
		switch (m_task) {
		case Task::PR:
		case Task::MAR:
			{
				double val = m_logz + std::log(orig.get_global_const());
				double prob = std::exp(val);

				out << "PR" << std::endl;
				out << std::fixed << std::setprecision(MERLIN_PRECISION)
					<< (m_logz + std::log(orig.get_global_const())) << " ("
					<< std::scientific << std::setprecision(MERLIN_PRECISION)
					<< std::exp(m_logz + std::log(orig.get_global_const()))
					<< ")" << std::endl;

				if (prob == 0.0) {
					out << "STATUS" << std::endl;
					out << "false: Inconsistent evidence or underflow" << std::endl;
				} else {
					out << "STATUS" << std::endl;
					out << "true: Consistent evidence" << std::endl;
				}

				out << "MAR" << std::endl;
				out << orig.nvar() - dummies.size();
				for (vindex i = 0; i < orig.nvar(); ++i) {
					if (dummies.find(i) != dummies.end()) {
						continue; // skip dummy variables from output
					}

					variable v = orig.var(i);
					try { // evidence variable
						size_t val = evidence.at(i);
						out << " " << v.states();
						for (size_t k = 0; k < v.states(); ++k) {
							out << " " << std::fixed
								<< std::setprecision(MERLIN_PRECISION)
								<< (k == val ? 1.0 : 0.0);
						}
					} catch(std::out_of_range& e) { // non-evidence variable
						vindex vx = old2new.at(i);
						variable VX = var(vx);
						out << " " << VX.states();
						for (size_t k = 0; k < VX.states(); ++k) {
							out << " " << std::fixed
								<< std::setprecision(MERLIN_PRECISION)
								<< belief(VX)[k];
						}
					}
				} // end for
				out << std::endl;

				if (!m_query.empty()) {
					variable_set scope = m_marginal.variables();
					out << "JOINT_MAR : [";
					for (size_t i = 0; i < m_query.size(); ++i) {
						out << m_query[i];
						if (i < m_query.size() - 1) {
							out << ",";
						}
					}
					out << "]" << std::endl;

					std::vector<size_t> dims(scope.size());
					size_t j = 0;
					for (variable_set::const_iterator si = scope.begin();
							si != scope.end(); ++si) {
						dims[j++] = (*si).states();
					}

					for (size_t index = 0; index < m_marginal.numel(); ++index) {
						std::vector<size_t> I(dims.size(), 0); // configuration of source variables
						size_t i = index;
						for (size_t v = 0; v < dims.size(); ++v) {
							I[v] = i % dims[v];
							i -= I[v];
							i /= dims[v];
						}

						for (size_t j = 0; j < I.size(); ++j) {
							out << I[j] << " ";
						}
						out << ": " << std::fixed
								<< std::setprecision(MERLIN_PRECISION)
								<< m_marginal[index]
								<< std::endl;
					}
				}

				break;
			}
		default:
			break;
		}
	} else {
		std::string err_msg("Unknown output format.");
		throw std::runtime_error(err_msg);
	}
}

} // namespace

ext: cpp

max_stars_repo_name: Nedriano/merlin

max_issues_repo_name: Nedriano/merlin

max_forks_repo_name: Nedriano/merlin


text:
#pragma once

#include <boost/exception/to_string.hpp>

namespace bunsan {
namespace process {

class context;

}  // namespace process
}  // namespace bunsan

namespace boost {
std::string to_string(const bunsan::process::context &context);
}  // namespace boost


text:
#include <QtGui>
#include <QtOpenGL>
#include <GL/glu.h>
#include <iostream>
#include <boost/foreach.hpp>
#include <cmath>
#include "TexShape.h"

#include "midasgl.h"

MidasGL :: MidasGL(vector< vector<Point> > p, vector<TexShape> tp, QTimer *timer, QWidget *parent)
	: QGLWidget(parent), polygons(p), texturedPolygons(tp) {
	trolltechGreen = QColor::fromCmykF(0.40,0.0,1.0,0.0);
	trolltechPurple = QColor::fromCmykF(0.39,0.39,0.0,0.0);
	mouseMoving = false;
	flightMode = true;
	this->timer = timer;
	frames = 0;
	initPos();
//	for (int i = 0; i < 10; i++) {
//
//	}
}

MidasGL::~MidasGL() {
	makeCurrent();
	glDeleteLists(polygonLists, polygons.size());
}

void MidasGL::initPos() {
	setXRotation(0);
	setYRotation(0);
	setZRotation(0);
	setEyeX(0);
	setEyeY(0);
	setEyeZ(350);
	setCenterX(0);
	setCenterY(0);
	setCenterZ(-100);
	yRotRoll = 0;
	zRotRoll = 0;
	xTrans = 0;
	yTrans = 0;
	zTrans = 0;
	upX = 0;
	upY = 1;
	upZ = 0;
	speed = 0.1;

}

QSize MidasGL::minimumSizeHint() const {
	return QSize(400,400);
}

QSize MidasGL::sizeHint() const {
	return QSize(400,400);
}

void MidasGL::setWireframe(int state) {
	switch (state) {
	case Qt::Unchecked:
		glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);
		break;
	case Qt::Checked:
		glPolygonMode(GL_FRONT_AND_BACK, GL_LINE);
		break;
	}
	updateGL();
}

void MidasGL::setSpeed(int speed) {
	this->speed = (double)speed / 10.0;
	emit speedChanged((int)(speed * 10));
}

void MidasGL::setXRotation(double angle) {
//	normalizeAngle(&angle);
	if (angle != xRot) {
		xRot = angle;
		emit xRotationChanged(angle);
//		setCenterY(eyeY + -200 * sin(xRot * (PI / 180.0)));
		updateGL();
	}
}
void MidasGL::setYRotation(double angle) {
//	normalizeAngle(&angle);
	if (angle != yRot) {
		yRot = angle;
		emit yRotationChanged(angle);
		setCenterX(eyeX + 200 * sin(yRot * (PI / 180.0)));
		setCenterZ(eyeZ + -200 * cos(yRot * (PI / 180.0)));
		//updateGL();
	}
}

void MidasGL::setZRotation(double angle) {
//	normalizeAngle(&angle);
	if (angle != zRot) {
		zRot = angle;
		emit zRotationChanged(angle);
		//updateGL();
	}
}
void MidasGL::setXRotation(int angle) {
	setXRotation((double)angle);
	emit xRotationChanged(angle);
}
void MidasGL::setYRotation(int angle) {
	setYRotation((double)angle);
	emit yRotationChanged(angle);
}
void MidasGL::setZRotation(int angle) {
	setZRotation((double)angle);
	emit yRotationChanged(angle);
}

void MidasGL::setXTrans(int xTrans) {
	this->xTrans = xTrans;
	//updateGL();
}
void MidasGL::setYTrans(int yTrans) {
	this->yTrans = yTrans;
	//updateGL();
}
void MidasGL::setZTrans(int zTrans) {
	this->zTrans = zTrans;
	if (!flightMode) updateGL();
}
void MidasGL::setEyeX(double eyeX) {
	this->eyeX = eyeX;
	emit xEyeChanged((int)eyeX);
	//updateGL();
}
void MidasGL::setEyeY(double eyeY) {
	this->eyeY = eyeY;
	emit yEyeChanged((int)eyeY);
	//updateGL();
}
void MidasGL::setEyeZ(double eyeZ) {
	this->eyeZ = eyeZ;
	emit zEyeChanged((int)eyeZ);
	if (!flightMode) updateGL();
}
void MidasGL::setEyeX(int eyeX) {
	setEyeX((double)eyeX);
}
void MidasGL::setEyeY(int eyeY) {
	setEyeY((double)eyeY);
}
void MidasGL::setEyeZ(int eyeZ) {
	setEyeZ((double)eyeZ);
}
void MidasGL::setCenterX(double centerX) {
	this->centerX = centerX;
	emit xCenterChanged((int)centerX);
	//updateGL();
}
void MidasGL::setCenterY(double centerY) {
	this->centerY = centerY;
	emit yCenterChanged((int)centerY);
	//updateGL();
}
void MidasGL::setCenterZ(double centerZ) {
	this->centerZ = centerZ;
	emit zCenterChanged((int)centerZ);
	//updateGL();
}
void MidasGL::setCenterX(int centerX) {
	setCenterX((double)centerX);
}
void MidasGL::setCenterY(int centerY) {
	setCenterY((double)centerY);
}
void MidasGL::setCenterZ(int centerZ) {
	setCenterZ((double)centerZ);
}
void MidasGL::setUpX(int upX) {
	this->upX = upX;
	//updateGL();
}
void MidasGL::setUpY(int upY) {
	this->upY = upY;
	//updateGL();
}
void MidasGL::setUpZ(int upZ) {
	this->upZ = upZ;
	//updateGL();
}

void MidasGL::flight() {
	this->timer->stop();
	if (flightMode) {
		elapsed = (elapsed + qobject_cast<QTimer*>(sender())->interval());
	//	cout << "elapsed: " << elapsed << endl;
	//	cout << "speed: " << speed << endl;
		realign();
		if (abs(zRot/16) > 1) { // rotate left/right a bit when rolling
			if (zRot > 1)
				setYRotation(yRot + log((double)zRot/32));
			else
				setYRotation(yRot - log((double)abs(zRot)/32));
		}
		if (abs(xRot/16) > 1) { // move the center a bit when mouse cause x rotation
			setCenterY(eyeY - xRot/2);
		}
		Vector velocity = Point(centerX, centerY, centerZ) - Point(eyeX, eyeY, eyeZ);
		velocity.normalize();
	//	cout << "direction: " << velocity << endl;
	//	cout << "velocity: " << velocity*speed << endl;
		double eyeXMove = eyeX + velocity.x()*speed;
		double eyeYMove = eyeY + velocity.y()*speed;
		double eyeZMove = eyeZ + velocity.z()*speed;
		setEyeX(eyeXMove);
		setEyeY(eyeYMove);
		setEyeZ(eyeZMove);
		setCenterX(centerX + velocity.x()*speed);
		setCenterY(centerY + velocity.y()*speed);
		setCenterZ(centerZ + velocity.z()*speed);
	}
	updateGL();
	frames++;
	this->timer->start();
}

void MidasGL::realign() {
	if (!mouseMoving && flightMode) {
		if (zRot != 0) {
			cout << "realigning zRot...\n";
			int zRotBack = (int)(.8 * (double)zRot);
			setZRotation(zRotBack);
			cout << "zRot: " << zRot << endl;
		}
		if (xRot != 0) {
//			cout << "realigning xRot...\n";
//			int xRotBack = (int)(.8 * (double)xRot);
//			setXRotation(xRotBack);
//			cout << "xRot: " << xRot << endl;
			setXRotation(0);
		}
		if (centerY != eyeY) {
			// if view is not leveled
			cout << "realigning centerY...\n";
			cout << "centerY: " << centerY << " eyeY: " << eyeY << endl;
			double yLevelBack = .2 * (centerY - eyeY);
			setCenterY(centerY - yLevelBack);
			if (abs(centerY - eyeY) <= 1)
				setCenterY(eyeY);
		}
	}
}

void MidasGL::initializeGL() {

	qglClearColor(trolltechPurple.dark());
	createPolygonLists();

//	const GLfloat red[4] = {0.6f, 0.f, 0.f, 1.f};
//	const GLfloat blue[4] = {0.f, 0.f, 0.4f, 1.f };
//	const GLfloat green[4] = {0.f, 0.5f, 0.f, 1.f };
	const GLfloat white[4] = {0.9f, 0.9f, 0.9f, 0.9f};
//	const GLfloat LightAmbient[] = { 0.5f, 0.5f, 0.5f, 1.0f };
	const GLfloat position[4] = { 0.f, 5.f, 20.f, 1.f };

	glLightModeli(GL_LIGHT_MODEL_TWO_SIDE, 0);
	glColorMaterial(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE);
	glShadeModel(GL_SMOOTH);
	glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);

	glLightfv( GL_LIGHT1, GL_AMBIENT, white);
	glLightfv( GL_LIGHT1, GL_DIFFUSE, white);
	glLightfv( GL_LIGHT1, GL_SPECULAR, white);
	glLightfv( GL_LIGHT1, GL_POSITION, position);
	glLightf(GL_LIGHT1, GL_CONSTANT_ATTENUATION, 2.f);
	glLightf(GL_LIGHT1, GL_LINEAR_ATTENUATION, 0.f);
	glLightf(GL_LIGHT1, GL_QUADRATIC_ATTENUATION, 0.f);

	glEnable(GL_DEPTH_TEST);
	glEnable(GL_LIGHTING);
	glEnable(GL_COLOR_MATERIAL);
	glEnable(GL_LIGHT1);
 	glEnable(GL_TEXTURE_2D);
}

void MidasGL::paintGL() {
	glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
	glLoadIdentity();

	if (flightMode) {
		flightView();
	} else {
		mouseView();
	}

//	const GLfloat position[4] = { eyeX, eyeY, eyeZ, 0.f };
//	glLightfv( GL_LIGHT1, GL_POSITION, position);

	glListBase(polygonLists);
	GLuint lists[polygons.size()];
	for (unsigned int i = 0; i < polygons.size(); i++) {
		lists[i] = i;
	}
	glCallLists((GLsizei)polygons.size(), GL_UNSIGNED_INT, lists);

//	glListBase(texturedPolygonLists);
//	unsigned int tsize = 0;
//	BOOST_FOREACH(TexShape ts, texturedPolygons) {
//		tsize += ts.shape.size();
//	}
//	GLuint tlists[tsize];
//	for (unsigned int i = 0; i < tsize; i++) {
//		tlists[i] = i;
//	}
//	glCallLists((GLsizei)tsize, GL_UNSIGNED_INT, tlists);

	int base = texturedPolygonLists;
	for (int i = 0; i < texturedPolygons.size(); i++) {
		glBindTexture(GL_TEXTURE_2D, textureLists[i]);
		glListBase(base);
		GLuint lists[texturedPolygons[i].shape.size()];
		for (unsigned int j = 0; j < texturedPolygons[i].shape.size(); j++)
			lists[j] = j;
		glCallLists((GLsizei)texturedPolygons[i].shape.size(), GL_UNSIGNED_INT, lists);
		base += texturedPolygons[i].shape.size();
	}

}

void MidasGL::mouseView() {
	glRotated(xRot/16.0, 1.0,0.0,0.0); // for mouse view
	glRotated(yRot/16.0, 0.0,1.0,0.0); //
	glRotated(zRot/16.0,0.0,0.0,1.0);  //
}

void MidasGL::flightView() {
	glTranslated(xTrans, yTrans, zTrans);
//	glRotated(xRot/16, 1.0,0.0,0.0); // for mouse view
//	glRotated(yRot, 0.0,1.0,0.0); //
	glRotated(zRot/16,0.0,0.0,1.0);  //
	gluLookAt(eyeX, eyeY, eyeZ, centerX, centerY, centerZ, upX, upY, upZ);

}

void MidasGL::resizeGL(int width, int height) {
	this->width = width;
	this->height = height;
	int side = qMin(width, height);
	glViewport((width-side)/2, (height-side)/2,side,side);
	glMatrixMode(GL_PROJECTION);
	glLoadIdentity();
	if (!flightMode) {
		glOrtho(-100.0, +100.0, -100.0, +100.0, -250.0, +250.0);
	} else {
		gluPerspective(45, width/height, 1, 550);
	}
	glMatrixMode(GL_MODELVIEW);
}

void MidasGL::mousePressEvent(QMouseEvent *event) {
	mouseMoving = true;
	lastPos = event->pos();
#if QT_VERSION < 0x040600
	if (event->buttons() & Qt::MiddleButton) {
#else
	if (event->buttons() & Qt::MidButton) {
#endif
		flightMode ^= true;
		cout << "Flight mode: " << flightMode << endl;
		if (flightMode) {
			initPos();
		}
		resizeGL(this->width, this->height);
		updateGL();
	}
}

void MidasGL::mouseReleaseEvent(QMouseEvent *event) {
	event->isAccepted();
	mouseMoving = false;
}

void MidasGL::wheelEvent(QWheelEvent *event) {
	int numDegrees = event->delta() / 8;
	int numSteps = numDegrees/15;
	double speedDelta = (double)numSteps / 10.0;
	speed += speedDelta;
	if (speed < 0)
		speed = 0;
	else if (speed > 5)
		speed = 5;
	emit speedChanged((int)(speed*10));
}

void MidasGL::mouseMoveEvent(QMouseEvent *event) {
	mouseMoving = true;
	int dx = event->x() - lastPos.x();
	int dy = event->y() - lastPos.y();

	if (event->buttons() & Qt::LeftButton) {
		if (flightMode) {
			double xRotation = xRot + dy;
			if (xRotation < -1440) xRotation = -1440;
			if (xRotation > 1440) xRotation = 1440;
			setXRotation(xRotation);
			//setYRotation(yRot + dx);

	//		int zRotate = zRot + 2 * dx;
	//		if (zRotate > 720) zRotate = 720;
	//		if (zRotate < -720) zRotate = -720;
	//		setZRotation(zRotate);
	//		zRotRoll += 2 * dx;
		} else {
			setXRotation(xRot + 8 * dy);
			setYRotation(yRot + 8 * dx);
		}
	}
	if (event->buttons() & Qt::RightButton) {
		if (flightMode) {
	//		int yRotate = yRot + dx;
	//		if (yRotate > 45) yRotate = 45;
	//		if (yRotate < -45) yRotate = -45;
	//		setYRotation(yRotate);
	//		if (yRot + dx <= 45 && yRot + dx >= -45)
	//			yRotRoll += dx;
			int zRotate = zRot + 2 * dx;
			if (zRotate > 720) zRotate = 720;
			if (zRotate < -720) zRotate = -720;
			setZRotation(zRotate);
			zRotRoll += 2 * dx;
		} else {
			setXRotation(xRot + 8 * dy);
			setZRotation(zRot + 8 * dx);
		}
	}
	cout << "xRot: " << xRot << endl;
	cout << "yRot: " << yRot << endl;
	cout << "zRot: " << zRot << endl;

	lastPos = event->pos();
}

void MidasGL::createPolygonLists() {
	polygonLists = glGenLists(polygons.size());
	for (unsigned int i = 0; i < polygons.size(); i++) {
		glNewList(polygonLists+i, GL_COMPILE);
		glBegin(GL_POLYGON);
		vector<Point> polygon = polygons[i];
		Vector vnorm;
		if (polygon.size() > 2) {
			Vector v1 = polygon.back() - polygon[0];
			Vector v2 = polygon[1] - polygon[0];
			vnorm = cross(v1, v2);
//			Vector vorigin = Point(0,0,0) - polygon[0];
//			if (dot(vnorm, vorigin) >= 0.0)
//				vnorm = -1*vnorm;
			vnorm.normalize();
			glNormal3f(vnorm.x(), vnorm.y(), vnorm.z());
		}
		BOOST_FOREACH(Point p, polygon) {
			GLubyte r = p.r(), g = p.g(), b = p.b(), a = p.a();
			glColor4ub(r,g,b,a);
			//glColor4f(p.r(), p.g(), p.b(), p.a());
			glVertex3f(p.x(), p.y(), p.z());
		}
		glEnd();
		glFlush();
		glEndList();
	}

	// Create texture handles.
	textureLists = new GLuint[texturedPolygons.size()];

	for (unsigned int i = 0; i < texturedPolygons.size(); i++) {
		// Load the image using Qt (can be done countless other ways too).
		QImage image(string("textures/").append(texturedPolygons[i].texturefile).c_str());
		QImage ogl = QGLWidget::convertToGLFormat(image);
		//cout << texturedPolygons[i].texturefile << endl;
		//cout << ogl.width() << " " << ogl.height() << endl;

		glGenTextures(1, &textureLists[i]);
		// Bind the texture.
		glBindTexture(GL_TEXTURE_2D, textureLists[i]);

		// Setup the texture parameters (important!).
		glTexEnvi(GL_TEXTURE_ENV, GL_TEXTURE_ENV_MODE, GL_REPLACE);
		glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_S, GL_CLAMP);
		glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_T, GL_CLAMP);
		glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_NEAREST);
		glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_NEAREST);
		glTexEnvf(GL_TEXTURE_ENV, GL_TEXTURE_ENV_MODE, GL_MODULATE);

		// Upload the texture to OpenGL.
		glTexImage2D(GL_TEXTURE_2D,
					 0,
					 GL_RGBA,
					 ogl.width(),
					 ogl.height(),
					 0,
					 GL_RGBA,
					 GL_UNSIGNED_BYTE,
					 ogl.bits());

		// Unbind the texture.
		glBindTexture(GL_TEXTURE_2D, 0);
	}
//	for (int a = 0; a < texturedPolygons.size(); a++) {
//		cout << textureLists[a] << endl;
//	}
//	exit(1);

	int tsize = 0;
	BOOST_FOREACH(TexShape ts, texturedPolygons) {
		tsize += ts.shape.size();
	}
	texturedPolygonLists = glGenLists(tsize);
	int listIndex = 0;
	for (unsigned int j = 0; j < texturedPolygons.size(); j++) {
		TexShape ts = texturedPolygons[j];
		GLint name;
		glGetIntegerv(GL_TEXTURE_BINDING_2D, &name);
		//cout << "current texture: " << name << endl;
		for(unsigned int i = 0; i < ts.shape.size(); i++) {
			glNewList(texturedPolygonLists+listIndex, GL_COMPILE);
			glBegin(GL_QUADS);

			vector<Point> polygon = ts.shape[i];
			vector<Point> texture = ts.texture[i];

			Vector vnorm;
			if (polygon.size() > 2) {
				Vector v1 = polygon.back() - polygon[0];
				Vector v2 = polygon[1] - polygon[0];
				vnorm = cross(v1, v2);
				vnorm.normalize();
				glNormal3f(vnorm.x(), vnorm.y(), vnorm.z());
			}
			for (unsigned int n = 0; n < polygon.size(); n++) {
				Point p = polygon[n];
				Point t = texture[n];
				GLubyte r = p.r(), g = p.g(), b = p.b(), a = p.a();
				GLfloat u = t.x(), v = t.y();
				glColor4ub(r,g,b,a);
				glTexCoord2f(u, v);
				glVertex3f(p.x(), p.y(), p.z());
			}
			glEnd();
			glFlush();
			glEndList();
			listIndex++;
		}
	}
	//cout << "tsize: " << tsize << endl;
	//cout << "listIndex: " << listIndex << endl;

}

void MidasGL::normalizeAngle(int *angle)
{
	while (*angle < 0)
		*angle += 360*16;
	while (*angle > 360*16)
		*angle -= 360*16;
}

max_stars_repo_name: eddy-chan/midas

max_issues_repo_name: eddy-chan/midas

max_forks_repo_name: eddy-chan/midas


text:
//////////////////////////////////////////////////////////////////////////////
//
// (C) Copyright Ion Gaztanaga 2005-2012. 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/interprocess for documentation.
//
//////////////////////////////////////////////////////////////////////////////

#ifndef BOOST_INTERPROCESS_DETAIL_ADAPTIVE_NODE_POOL_HPP
#define BOOST_INTERPROCESS_DETAIL_ADAPTIVE_NODE_POOL_HPP

#ifndef BOOST_CONFIG_HPP
#  include <boost/config.hpp>
#endif
#
#if defined(BOOST_HAS_PRAGMA_ONCE)
#  pragma once
#endif

#include <boost/interprocess/detail/config_begin.hpp>
#include <boost/interprocess/detail/workaround.hpp>
#include <boost/interprocess/detail/utilities.hpp>
#include <boost/interprocess/detail/math_functions.hpp>
#include <boost/intrusive/set.hpp>
#include <boost/intrusive/slist.hpp>
#include <boost/interprocess/detail/type_traits.hpp>
#include <boost/interprocess/mem_algo/detail/mem_algo_common.hpp>
#include <boost/interprocess/allocators/detail/node_tools.hpp>
#include <boost/interprocess/allocators/detail/allocator_common.hpp>
#include <cstddef>
#include <boost/config/no_tr1/cmath.hpp>
#include <boost/container/detail/adaptive_node_pool_impl.hpp>
#include <boost/assert.hpp>

//!\file
//!Describes the real adaptive pool shared by many Interprocess pool allocators

namespace boost {
namespace interprocess {
namespace ipcdetail {

template< class SegmentManager
        , std::size_t NodeSize
        , std::size_t NodesPerBlock
        , std::size_t MaxFreeBlocks
        , unsigned char OverheadPercent
        >
class private_adaptive_node_pool
   :  public boost::container::container_detail::private_adaptive_node_pool_impl
         < typename SegmentManager::segment_manager_base_type
         , ::boost::container::adaptive_pool_flag::size_ordered |
           ::boost::container::adaptive_pool_flag::address_ordered
         >
{
   typedef boost::container::container_detail::private_adaptive_node_pool_impl
      < typename SegmentManager::segment_manager_base_type
      , ::boost::container::adaptive_pool_flag::size_ordered |
        ::boost::container::adaptive_pool_flag::address_ordered
      > base_t;
   //Non-copyable
   private_adaptive_node_pool();
   private_adaptive_node_pool(const private_adaptive_node_pool &);
   private_adaptive_node_pool &operator=(const private_adaptive_node_pool &);

   public:
   typedef SegmentManager              segment_manager;
   typedef typename base_t::size_type  size_type;

   static const size_type nodes_per_block = NodesPerBlock;

   //Deprecated, use node_per_block
   static const size_type nodes_per_chunk = NodesPerBlock;

   //!Constructor from a segment manager. Never throws
   private_adaptive_node_pool(segment_manager *segment_mngr)
      :  base_t(segment_mngr, NodeSize, NodesPerBlock, MaxFreeBlocks, OverheadPercent)
   {}

   //!Returns the segment manager. Never throws
   segment_manager* get_segment_manager() const
   {  return static_cast<segment_manager*>(base_t::get_segment_manager_base()); }
};

//!Pooled shared memory allocator using adaptive pool. Includes
//!a reference count but the class does not delete itself, this is
//!responsibility of user classes. Node size (NodeSize) and the number of
//!nodes allocated per block (NodesPerBlock) are known at compile time
template< class SegmentManager
        , std::size_t NodeSize
        , std::size_t NodesPerBlock
        , std::size_t MaxFreeBlocks
        , unsigned char OverheadPercent
        >
class shared_adaptive_node_pool
   :  public ipcdetail::shared_pool_impl
      < private_adaptive_node_pool
         <SegmentManager, NodeSize, NodesPerBlock, MaxFreeBlocks, OverheadPercent>
      >
{
   typedef ipcdetail::shared_pool_impl
      < private_adaptive_node_pool
         <SegmentManager, NodeSize, NodesPerBlock, MaxFreeBlocks, OverheadPercent>
      > base_t;
   public:
   shared_adaptive_node_pool(SegmentManager *segment_mgnr)
      : base_t(segment_mgnr)
   {}
};

}  //namespace ipcdetail {
}  //namespace interprocess {
}  //namespace boost {

#include <boost/interprocess/detail/config_end.hpp>

#endif   //#ifndef BOOST_INTERPROCESS_DETAIL_ADAPTIVE_NODE_POOL_HPP

max_stars_repo_name: multi-os-engine/cinder-natj-binding


text:
// Copyright David Abrahams 2003. Use, modification and distribution is
// subject to the Boost Software License, Version 1.0. (See accompanying
// file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
#ifndef IS_LVALUE_ITERATOR_DWA2003112_HPP
# define IS_LVALUE_ITERATOR_DWA2003112_HPP

#include <boost/iterator.hpp>

#include <boost/detail/workaround.hpp>
#include <boost/detail/iterator.hpp>

#include <boost/iterator/detail/any_conversion_eater.hpp>

// should be the last #includes
#include <boost/type_traits/detail/bool_trait_def.hpp>
#include <boost/iterator/detail/config_def.hpp>

#ifndef BOOST_NO_IS_CONVERTIBLE

namespace boost {

namespace detail
{
#ifndef BOOST_NO_LVALUE_RETURN_DETECTION
  // Calling lvalue_preserver( <expression>, 0 ) returns a reference
  // to the expression's result if <expression> is an lvalue, or
  // not_an_lvalue() otherwise.
  struct not_an_lvalue {};

  template <class T>
  T& lvalue_preserver(T&, int);

  template <class U>
  not_an_lvalue lvalue_preserver(U const&, ...);

# define BOOST_LVALUE_PRESERVER(expr) detail::lvalue_preserver(expr,0)

#else

# define BOOST_LVALUE_PRESERVER(expr) expr

#endif

  // Guts of is_lvalue_iterator.  Value is the iterator's value_type
  // and the result is computed in the nested rebind template.
  template <class Value>
  struct is_lvalue_iterator_impl
  {
      // Eat implicit conversions so we don't report true for things
      // convertible to Value const&
      struct conversion_eater
      {
          conversion_eater(Value&);
      };

      static char tester(conversion_eater, int);
      static char (& tester(any_conversion_eater, ...) )[2];

      template <class It>
      struct rebind
      {
          static It& x;

          BOOST_STATIC_CONSTANT(
              bool
            , value = (
                sizeof(
                    is_lvalue_iterator_impl<Value>::tester(
                        BOOST_LVALUE_PRESERVER(*x), 0
                    )
                ) == 1
            )
          );
      };
  };

#undef BOOST_LVALUE_PRESERVER

  //
  // void specializations to handle std input and output iterators
  //
  template <>
  struct is_lvalue_iterator_impl<void>
  {
      template <class It>
      struct rebind : boost::mpl::false_
      {};
  };

#ifndef BOOST_NO_CV_VOID_SPECIALIZATIONS
  template <>
  struct is_lvalue_iterator_impl<const void>
  {
      template <class It>
      struct rebind : boost::mpl::false_
      {};
  };

  template <>
  struct is_lvalue_iterator_impl<volatile void>
  {
      template <class It>
      struct rebind : boost::mpl::false_
      {};
  };

  template <>
  struct is_lvalue_iterator_impl<const volatile void>
  {
      template <class It>
      struct rebind : boost::mpl::false_
      {};
  };
#endif

  //
  // This level of dispatching is required for Borland.  We might save
  // an instantiation by removing it for others.
  //
  template <class It>
  struct is_readable_lvalue_iterator_impl
    : is_lvalue_iterator_impl<
          BOOST_DEDUCED_TYPENAME boost::detail::iterator_traits<It>::value_type const
      >::template rebind<It>
  {};

  template <class It>
  struct is_non_const_lvalue_iterator_impl
    : is_lvalue_iterator_impl<
          BOOST_DEDUCED_TYPENAME boost::detail::iterator_traits<It>::value_type
      >::template rebind<It>
  {};
} // namespace detail

// Define the trait with full mpl lambda capability and various broken
// compiler workarounds
BOOST_TT_AUX_BOOL_TRAIT_DEF1(
    is_lvalue_iterator,T,::boost::detail::is_readable_lvalue_iterator_impl<T>::value)

BOOST_TT_AUX_BOOL_TRAIT_DEF1(
    is_non_const_lvalue_iterator,T,::boost::detail::is_non_const_lvalue_iterator_impl<T>::value)

} // namespace boost

#endif

#include <boost/iterator/detail/config_undef.hpp>
#include <boost/type_traits/detail/bool_trait_undef.hpp>

#endif // IS_LVALUE_ITERATOR_DWA2003112_HPP

max_stars_repo_name: randolphwong/mcsema

max_issues_repo_name: randolphwong/mcsema

max_forks_repo_name: randolphwong/mcsema


text:
//   GAMBIT: Global and Modular BSM Inference Tool
//   *********************************************
///  \file
///
///  Type definition header for module FlavBit.
///
///  Compile-time registration of type definitions
///  required for the rest of the code to
///  communicate with FlavBit.
///
///  Add to this if you want to define a new type
///  for the functions in FlavBit to return, but
///  you don't expect that type to be needed by
///  any other modules.
///
///  *********************************************
///
///  Authors (add name and date if you modify):
///
///  \author Marcin Chrzaszcz
///          (mchrzasz@cern.ch)
///  \date 2015 July
///
///  \author Tomas Gonzalo
///          (t.e.gonzalo@fys.uio.no)
///  \date 2017 July
///
///  *********************************************

#ifndef __FlavBit_types_hpp__
#define __FlavBit_types_hpp__

#include <string>
#include <vector>

#include <boost/numeric/ublas/matrix.hpp>

namespace Gambit
{

  namespace FlavBit
  {

    namespace ublas = boost::numeric::ublas;

    /// Simple structure for holding a correlation value and name of the correlated observable
    struct Correlation
    {
      double corr_val;
      std::string corr_name;
    };

    /// Representation of a single entry in the FlavBit YAML database
    struct Measurement
    {
      bool is_limit;
      bool exp_one_sided;
      double exp_value;
      double exp_stat_error;
      double exp_sys_error;
      double exp_error;
      double exp_90_CL;
      double exp_95_CL;
      double th_error;
      std::string name;
      std::string exp_source;
      std::string th_error_source;
      std::string th_error_type;
      std::vector<Correlation> corr;
    };

    /// Structure for holding predicted and observed values of multiple observables,
    /// and experimental and theoretical covariance matrices for their uncertainties.
    struct predictions_measurements_covariances
    {
      std::string LL_name;
      boost::numeric::ublas::matrix<double> value_exp;
      boost::numeric::ublas::matrix<double> cov_exp;
      boost::numeric::ublas::matrix<double> value_th;
      boost::numeric::ublas::matrix<double> cov_th;
      std::vector<double> diff;
      int dim;
    };

  }
}

#endif


text:
#pragma once

#include <memory>
#include <iostream>

#include <boost/asio.hpp>
#include <boost/bind.hpp>

#include "json.hpp"

using namespace boost;

class Connection;
class ServicePool;

class WebService {
public:
	std::map<std::string, std::string(*)(JsonObject* json)> routes;
	std::map<std::string, std::vector<std::pair<std::string, JsonType>>> required_fields;

	WebService(int port, std::string service_name);
	void route(std::string path, std::string(*func)(JsonObject* json), std::vector<std::pair<std::string, JsonType>> requires = {});
	void listen();
	std::string GetName();
private:
	int port;
	std::string service_name;
	ServicePool* service_pool;
	asio::ip::tcp::acceptor serv_acceptor;

	void connection_acceptor(const system::error_code& ec, Connection* new_connection);
};

max_stars_repo_name: bwackwat/raw-json-webservices

max_issues_repo_name: bwackwat/cpp-json-webservices

max_forks_repo_name: bwackwat/cpp-json-webservices


text:
#ifndef boost_xml_backends_xerces_reader_hpp_
#define boost_xml_backends_xerces_reader_hpp_

#include <boost/xml/converter.hpp>
#include <xercesc/framework/XMLPScanToken.hpp>
#include <xercesc/parsers/SAXParser.hpp>
#include <xercesc/sax/HandlerBase.hpp>
#include <stdexcept>
#include <iostream>

namespace boost
{
namespace xml
{
namespace xerces
{
namespace reader
{
enum token_type
{
  none,
  element,
  attribute,
  text,
  cdata,
  entity_reference,
  type_entity,
  processing_instruction,
  comment,
  document,
  document_type,
  document_fragment,
  notation,
  whitespace,
  significant_whitespace,
  end_element,
  end_entity,
  xml_declaration,
};

XERCES_CPP_NAMESPACE_USE

template <typename S>
class handler : public HandlerBase
{
public :
  handler() : type_(none), depth_(0) {}
  ~handler() {}

  int depth() const { return depth_;}
  token_type type() const { return type_;}
  XMLCh const *name() const { return name_;}
  XMLCh const *value() const { return value_;}

  void resetDocType()
  {
    type_ = document_type;
  }
  void notationDecl(XMLCh const *const name, XMLCh const *const publicId, XMLCh const *const systemId)
  {
    type_ = document_type;
    name_ = name;
  }
  void startElement(XMLCh const *const name, AttributeList &attributes)
  {
    type_ = element;
    ++depth_;
    name_ = name;
  }
  void endElement(XMLCh const *const name) { --depth_;}
  void characters(XMLCh const *const chars, XMLSize_t length)
  {
    type_ = text;
    name_ = converter<S>::in("#text");
    value_ = chars;
    length_ = length;
  }
  void ignorableWhitespace(XMLCh const *const chars, XMLSize_t length)
  {
    type_ = whitespace;
    value_ = chars;
    length_ = length;
  }
  void processingInstruction(XMLCh const *target, XMLCh const *data)
  {
    type_ = processing_instruction;
    name_ = target;
    value_ = data;
  }
  void resetDocument()
  {
  }

#if 0
  void warning(SAXParseException const &e)
  {
    std::cerr << "\nWarning at file " << converter<S>::out(e.getSystemId())
	      << ", line " << e.getLineNumber()
	      << ", char " << e.getColumnNumber()
	      << "\n  Message: " << converter<S>::out(e.getMessage()) << std::endl;
  }
  void error(SAXParseException const &e)
  {
    std::cerr << "\nError at file " << converter<S>::out(e.getSystemId())
	      << ", line " << e.getLineNumber()
	      << ", char " << e.getColumnNumber()
	      << "\n  Message: " << converter<S>::out(e.getMessage()) << std::endl;
  }
  void fatalError(SAXParseException const &e)
  {
    std::cerr << "\nFatal Error at file " << converter<S>::out(e.getSystemId())
	      << ", line " << e.getLineNumber()
	      << ", char " << e.getColumnNumber()
	      << "\n  Message: " << converter<S>::out(e.getMessage()) << std::endl;
  }
#endif

private:
  token_type type_;
  int depth_;
  XMLCh const *name_;
  XMLCh const *value_;
  XMLSize_t length_;
};

template <typename S> class parser;

template <typename S>
class token
{
  friend class parser<S>;
public:
  size_t depth() const { return handler_.depth();}
  S name() const { return converter<S>::out(handler_.name());}
  S value() const { return converter<S>::out(handler_.value());}
  bool is_empty_element() const { return handler_.is_empty_element();}
  bool has_value() const { return handler_.value() != 0;}
  token_type type() const { return handler_.type();}
protected:
  token(handler<S> &h) : handler_(h) {}

  handler<S> &handler_;
};

template <typename S>
class parser
{
public:
  parser(std::string const &filename)
    : first_(true)
  {
    parser_.setDTDHandler(&handler_);
    parser_.setDocumentHandler(&handler_);
    if (!parser_.parseFirst(filename.c_str(), token_))
      throw std::runtime_error("unable to open file");
  }
  ~parser() {}
  bool next()
  {
    if (first_)
    {
      first_ = false;
      return true;
    }
    else
      return parser_.parseNext(token_);
  }
  reader::token<S> token()
  {
    reader::token<S> t(handler_);
    return t;
  }
private:
  handler<S> handler_;
  SAXParser parser_;
  XMLPScanToken token_;
  bool first_;
};
}
}
}
}

#endif


text:
/**
 * @file k_means_clustering.hpp
 * @brief
 * @author Piotr Smulewicz, Piotr Wygocki
 * @version 1.0
 * @date 2014-06-25
 */
#ifndef PAAL_K_MEANS_CLUSTERING_HPP
#define PAAL_K_MEANS_CLUSTERING_HPP

#include "paal/clustering/k_means_clustering_engine.hpp"
#include "paal/utils/irange.hpp"

#include <boost/range/algorithm_ext/iota.hpp>
#include <boost/range/empty.hpp>

#include <vector>

namespace paal {

/**
 * @brief  return centroid that minimize within-cluster sum of squares
 */
template <typename Cluster, typename OutputIterator>
void centroid_minimalize_w_c_s_s(Cluster && cluster, OutputIterator out) {
    assert(!boost::empty(cluster));
    using point_t = range_to_elem_t<Cluster>;
    using coordinate_t = range_to_elem_t<point_t>;

    auto dim = boost::size(*std::begin(cluster));
    for(auto idx : irange(dim)) {
        coordinate_t res{};
        for (auto && point : cluster) {
            res += std::begin(point)[idx];
        }
        *out = res / boost::size(cluster);
        ++out;
   }
}

/**
 * @brief centroid minimize within cluster sum of squares
 * @param clusters
 * @param out
 * @tparam Clusters
 * @tparam OutputIterator
 */
template <typename Clusters, typename OutputIterator>
void centroids_minimalize_w_c_s_s(Clusters && clusters, OutputIterator out) {
    assert(!boost::empty(clusters));
    assert(!boost::empty(*std::begin(clusters)));

    using cluster_t = range_to_elem_t<Clusters>;
    using point_t = range_to_elem_t<cluster_t>;
    using coordinate_t = range_to_elem_t<point_t>;

    auto size = boost::size(*std::begin(*begin(clusters)));
    for (auto && cluster : clusters) {
        std::vector<coordinate_t> point(size);
        centroid_minimalize_w_c_s_s(cluster, point.begin());
        *out = point;
        ++out;
    }
}

/**
 * @brief this is solve k_means_clustering problem
 * and return vector of cluster
 * example:
 *  \snippet k_means_clustering_example.cpp K Means Clustering Example
 *
 * complete example is k_means_clustering_example.cpp
 * @param points
 * @param centers
 * @param result pairs of point and id of cluster
 * (number form 0,1,2 ...,number_of_cluster-1)
 * @param visitor
 * @tparam Points
 * @tparam OutputIterator
 * @tparam CoordinateType
 * @tparam Visitor
 */
template <typename Points, class Centers, class OutputIterator,
          class Visitor = k_means_visitor>
auto k_means(Points &&points, Centers &&centers, OutputIterator result,
             Visitor visitor = Visitor{}) {
    using point_t = range_to_elem_t<Points>;
    using center_t = range_to_elem_t<Centers>;

    center_t center{ *std::begin(centers) };

    return k_means(
        points, centers,
        [&](std::vector<point_t> const & points)->center_t const & {
            centroid_minimalize_w_c_s_s(points, std::begin(center));
            return center;
        },
        [&](point_t const &point) { return closest_to(point, centers); },
        result, utils::equal_to{}, visitor);
}

} //!paal

#endif /* PAAL_K_MEANS_CLUSTERING_HPP */

max_stars_repo_name: Kommeren/AA

max_issues_repo_name: Kommeren/AA

max_forks_repo_name: Kommeren/AA