text:
// Author(s): Wieger Wesselink
// Copyright: see the accompanying file COPYING or copy at
// https://github.com/mCRL2org/mCRL2/blob/master/COPYING
//
// 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)
//
/// \file table_test.cpp
/// \brief Add your file description here.

#include <iostream>
#include <string>
#include <boost/test/minimal.hpp>

#include <map>
#include "mcrl2/atermpp/aterm_io.h"
#include "mcrl2/atermpp/indexed_set.h"

using namespace std;
using namespace atermpp;

void test_indexed_set()
{
  indexed_set<atermpp::aterm> t(100, 75);

  std::pair<long, bool> p;
  p = t.put(read_term_from_string("a"));
  BOOST_CHECK(t.size() == 1);
  p = t.put(read_term_from_string("b"));
  BOOST_CHECK(t.size() == 2);

  {
    const indexed_set<atermpp::aterm>& t1 = t;
  }
  indexed_set<atermpp::aterm> t2 = t;

  BOOST_CHECK(t.index(read_term_from_string("a")) == 0);
  BOOST_CHECK(t.index(read_term_from_string("b")) == 1);

  atermpp::aterm a = t.get(0);
  BOOST_CHECK(a == read_term_from_string("a"));

  atermpp::aterm b = t.get(1);
  BOOST_CHECK(b == read_term_from_string("b"));
  t.erase(a);
  BOOST_CHECK(t.size() == 1);

  t.erase(read_term_from_string("x"));
  BOOST_CHECK(t.size() == 1);

  p = t.put(read_term_from_string("c"));
  BOOST_CHECK(t.size() == 2);
  BOOST_CHECK(p.first == 0);

  t.clear();
  BOOST_CHECK(t.size() == 0);

  std::map<int, indexed_set<atermpp::aterm> > x;
  x[2] = t;
}

int test_main(int argc, char* argv[])
{
  test_indexed_set();
  return 0;
}

ext: cpp


text:
#include "network/sock_server.h"
#include "network/future.hpp"
#include <boost/log/trivial.hpp>
#include <mutex>
using std::size_t;

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

sock_server::sock_server(const config &conf, net_service *service)
    :
    port_(conf.this_node_config().port()),
    endpoint_(tcp::v4(), conf.this_node_config().port()),

    service_(service),
    stopped_(false),
    conf_(conf) {}

sock_server::~sock_server() {}

void sock_server::async_accept_connection() {
  ptr<tcp::socket> socket(
      new tcp::socket(service_->get_service(SERVICE_HANDLE)));
  ptr<connection> ctx = cs_new<connection>(socket, service_->get_handler(), false);
  acceptor_->async_accept(*socket, [this, ctx, socket](boost_ec ec) {
    if (stopped_.load()) {
      BOOST_LOG_TRIVIAL(info) << id_2_name(conf_.node_id()) << " stop accept connection";
      return;
    }
    if (not ec) {
      boost::asio::ip::tcp::endpoint ep = socket->remote_endpoint();
      std::string ip = ep.address().to_string() + ":" + std::to_string(ep.port());
      client_conn_mutex_.lock();
      incomming_conn_.insert(std::make_pair(ip, ctx));
      client_conn_mutex_.unlock();
      ctx->connected();
      async_accept_connection();
      async_accept_new_connection_done(ctx);
    } else {
      BOOST_LOG_TRIVIAL(info) << id_2_name(conf_.node_id()) << " async accept error ..";
      ctx->process_error(berror(ec));
    }
  });
}

void sock_server::async_accept_new_connection_done(
    ptr<connection> ctx) {
  ctx->async_read();
}

bool sock_server::start() {
  service_->start();

  acceptor_.reset(new boost::asio::ip::tcp::acceptor(service_->get_service(SERVICE_HANDLE),
                                                     endpoint_));
  socket_.reset(new boost::asio::ip::tcp::socket(service_->get_service(SERVICE_HANDLE)));

  async_accept_connection();
  BOOST_LOG_TRIVIAL(info) << id_2_name(conf_.node_id()) << " server listen on port [" << port_ << "].";
  BOOST_LOG_TRIVIAL(info) << id_2_name(conf_.node_id()) << " begin accept connection ...";
  return true;
}

void sock_server::stop() {
  if (stopped_.load()) {
    return;
  }
  client_conn_mutex_.lock();
  for (auto kv: incomming_conn_) {
    ptr<connection> c = kv.second;
    if (c) {
      c->close();
    }
  }
  incomming_conn_.clear();
  client_conn_mutex_.unlock();
  stopped_.store(true);
  if (service_) {
    service_->stop();
  }
  if (acceptor_) {
    acceptor_->close();
  }
  if (socket_) {
    socket_->close();
  }
}

void sock_server::join() {
  if (service_) {
    service_->join();
  }
}

max_stars_repo_name: ybbhwxfj/block-db

max_issues_repo_name: ybbhwxfj/block-db

max_forks_repo_name: ybbhwxfj/block-db


text:
/**
 * BCL to FASTQ file converter
 * Copyright (c) 2007-2017 Illumina, Inc.
 *
 * This software is covered by the accompanying EULA
 * and certain third party copyright/licenses, and any user of this
 * source file is bound by the terms therein.
 *
 * \file GzipCompressor.hh
 *
 * \brief Declaration of Gzip Compressor.
 *
 * \author Roman Petrovski
 * \author Aaron Day
 */

#ifndef BCL2FASTQ_IO_GZIP_COMPRESSOR_HH
#define BCL2FASTQ_IO_GZIP_COMPRESSOR_HH

#include <vector>
#include <boost/noncopyable.hpp>
#include <boost/iostreams/filtering_stream.hpp>
#include <boost/iostreams/filter/gzip.hpp>

namespace bcl2fastq
{
namespace io
{

namespace bios=boost::iostreams;

class GzipCompressor : private boost::noncopyable
{
public:
    GzipCompressor(std::vector<char>&       outputBuffer,
                   bool                     useBgzf,
                   const bios::gzip_params& gzipParams = bios::gzip::default_compression);

    ~GzipCompressor();

    std::streamsize write(const char* s, std::streamsize n);

    void flush();

    static const int SUBFIELD_LENGTH = 2;  // this needs to be 2 to comply with the BGZF specs

    static const unsigned MAX_BGZF_BLOCK_SIZE = 0xFFFF;

private:
#pragma pack(push, 1)
    struct Header
    {
        struct XFIELD
        {
            XFIELD(size_t compressedSize);

            unsigned char XLEN[2];

            unsigned char SI1;
            unsigned char SI2;
            unsigned char SLEN[2];
            unsigned char BSIZE[SUBFIELD_LENGTH];

            unsigned getXLEN() const {return unsigned(XLEN[0]) + XLEN[1] * 256;}
            unsigned getBSIZE() const {return unsigned(BSIZE[0]) + BSIZE[1] * 256;}

        };

        unsigned char ID1;
        unsigned char ID2;
        unsigned char CM;
        unsigned char FLG;
        unsigned char MTIME[4];
        unsigned char XFL;
        unsigned char OS;

        XFIELD xfield;
    };
#pragma pack(pop)

public:
    static size_t getHeaderSize() { return sizeof(Header); }

    static unsigned getBlockSize(char* buffer, bool failQuietly = false);

private:
    void prepareHeader();

    void rewriteHeader();

    // unverified number of bytes the compressor would add to data with compression level 0
    static const size_t gzip_junk = 41;
    // bgzf cannot handle blocks over 0xFFFF bytes long. assume 1/1 compression ratio for input data
    static const unsigned bgzf_buffer_size_ = MAX_BGZF_BLOCK_SIZE;
    static const unsigned short max_uncompressed_per_block_ = bgzf_buffer_size_ - gzip_junk;

    std::vector<char>&                            outputBuffer_;
    std::vector<char>                             internalBuffer_;
    bios::back_insert_device< std::vector<char> > sink_;
    boost::iostreams::gzip_compressor             compressor_;

    // as it is difficult to predict the size of compressed data, the buffering is based on
    // the amount of uncompressed data consumed. The assumption is that compressed data
    // will always require no more space than the uncompressed.
    size_t uncompressed_in_;

    size_t currentBlockBegin_;

    bool useBgzf_;
};

} // namespace io
} // namespace bcl2fastq

#endif // BCL2FASTQ_IO_GZIP_COMPRESSOR_HH

ext: hh


text:
// Copyright (c) 2019 Graphcore Ltd. All rights reserved.
#define BOOST_TEST_MODULE PipelineTrainingTest0

#include <../random_util.hpp>
#include <boost/test/unit_test.hpp>
#include <filereader.hpp>
#include <popart/builder.hpp>
#include <popart/dataflow.hpp>
#include <popart/devicemanager.hpp>
#include <popart/inputshapeinfo.hpp>
#include <popart/ndarraywrapper.hpp>
#include <popart/op/identity.hpp>
#include <popart/op/l1.hpp>
#include <popart/session.hpp>
#include <popart/sgd.hpp>
#include <popart/tensorinfo.hpp>
#include <popart/tensornames.hpp>
#include <popart/testdevice.hpp>

#include <algorithm>
#include <map>
#include <tuple>
#include <vector>

//
// In this depp linear model, the continuous weight update is equivalent to
// exact SGD. This enables us to test the final weights numerically.
//
BOOST_AUTO_TEST_CASE(ContinuousEquivalentTest0) {

  using namespace popart;
  int seed = 1011;
  DefaultRandomEngine eng(seed);
  UniformRealDistribution<float> fdis(0.f, 1.f);

  auto builder     = Builder::create();
  auto aiOnnx      = builder->aiOnnxOpset9();
  auto aiGraphcore = builder->aiGraphcoreOpset1();

  // accumulationFactor = replicationFactor = 1 in this example,
  // so batchSize = microBatchSize

  // samples in batch  (samples per weight update)
  int64_t batchSize = 4;

  int64_t batchesPerStep = 100;

  // samples in a step (samples processed with each call to run(...))
  int64_t samplesPerStep = batchesPerStep * batchSize;

  // the number of weight update steps is batchesPerStep = 100
  // the number of samples in a step is batchesPerStep * batchSize = 400

  int64_t sampleHeight = 8;
  std::vector<int64_t> sampleShape{sampleHeight, sampleHeight};

  // in this test, weights are of the same shape as data samples
  std::vector<int64_t> weightShape = sampleShape;

  // shape of a batch
  std::vector<int64_t> batchShape{batchSize, sampleHeight, sampleHeight};
  std::vector<int64_t> stepDataShape{
      batchesPerStep, batchSize, sampleHeight, sampleHeight};

  TensorInfo sampleInfo{"FLOAT", sampleShape};
  TensorInfo weightInfo = sampleInfo;
  TensorInfo batchInfo{"FLOAT", batchShape};
  TensorInfo stepDataInfo{"FLOAT", stepDataShape};

  int64_t sampleElms{sampleHeight * sampleHeight};
  int64_t batchElms    = sampleElms * batchSize;
  int64_t stepDataElms = batchElms * batchesPerStep;

  auto input1 = builder->addInputTensor(batchInfo, "1tupni");

  // The model : A series of additions of a Weight, and occassionally scaling
  // by 1.0 = (-1*-1).
  //
  // All weights are initialized as zero, and all inputs will be large (positive
  // or negative). This means that the l1 loss gradient (+lambda or -lambda)
  // will depend exclusively on the input (we use a small enough learning rate
  // so that this fact remains true during our training steps) This means that
  // the gradient applied to weights depends only on the sample for which it is
  // computed, which ensures the equivalence between continous and gradient
  // accumuated pipelining
  //
  std::vector<float> w0Vals(sampleElms, 0.0f);
  ConstVoidData w0Data = {w0Vals.data(), sampleInfo};
  auto w0              = builder->addInitializedInputTensor(w0Data);
  auto act0            = aiOnnx.add({w0, input1}, "act0");
  act0                 = aiGraphcore.scale({act0}, -1.0);
  act0                 = aiGraphcore.scale({act0}, -1.0);

  std::vector<float> w1Vals(sampleElms, 0.0f);
  ConstVoidData w1Data = {w1Vals.data(), sampleInfo};
  auto w1              = builder->addInitializedInputTensor(w1Data);
  auto act1            = aiOnnx.add({w1, act0}, "act1");

  std::vector<float> w2Vals(sampleElms, 0.0f);
  ConstVoidData w2Data = {w2Vals.data(), sampleInfo};
  auto w2              = builder->addInitializedInputTensor(w2Data);
  auto act2            = aiOnnx.add({w2, act1}, "act2");
  act2                 = aiGraphcore.scale({act2}, -1.0);
  act2                 = aiGraphcore.scale({act2}, -1.0);

  std::vector<float> w3Vals(sampleElms, 0.0f);
  ConstVoidData w3Data = {w3Vals.data(), sampleInfo};
  auto w3              = builder->addInitializedInputTensor(w3Data);
  auto act3            = aiOnnx.add({w3, act2}, "act3");

  std::vector<float> w4Vals(sampleElms, 0.0f);
  ConstVoidData w4Data = {w4Vals.data(), sampleInfo};
  auto w4              = builder->addInitializedInputTensor(w4Data);
  auto act4            = aiOnnx.add({w4, act3}, "act4");
  act4                 = aiGraphcore.scale({act4}, -1.0);
  act4                 = aiGraphcore.scale({act4}, -1.0);

  std::vector<float> w5Vals(sampleElms, 0.0f);
  ConstVoidData w5Data = {w5Vals.data(), sampleInfo};
  auto w5              = builder->addInitializedInputTensor(w5Data);
  auto act5            = aiOnnx.add({w5, act4}, "act5");
  float lambda         = 0.1;
  auto l1 =
      builder->aiGraphcoreOpset1().l1loss({act5}, lambda, ReductionType::Sum);

  auto proto = builder->getModelProto();

  // Setting anchors as act5
  auto dataFlow = DataFlow(batchesPerStep, {{act5, AnchorReturnType("All")}});

  // shard over 3 IPUs, and enable pipelining
  SessionOptions userOptions;
  userOptions.virtualGraphMode = VirtualGraphMode::Auto;
  userOptions.enablePipelining = true;
  std::map<std::string, std::string> deviceOpts{{"numIPUs", "3"}};

  float learnRate = 0.01;
  auto optimizer  = ConstSGD(learnRate);

  auto device = createTestDevice(TEST_TARGET, 3);

  auto session = popart::TrainingSession::createFromOnnxModel(
      proto,
      dataFlow,
      l1,
      optimizer,
      device,
      InputShapeInfo(),
      userOptions,
      popart::Patterns(PatternsLevel::Default));

  session->prepareDevice();

  // The 400 samples.
  std::vector<float> v_input_0(stepDataElms);

  // cumulative samples (accumulated over multiple steps).
  std::vector<float> v_sample_sum_0(weightInfo.nelms(), 0.0f);

  // Preparing the anchors.
  std::vector<float> v_act5_out(stepDataElms);
  popart::NDArrayWrapper<float> act5_wrapper(v_act5_out.data(), stepDataShape);

  std::map<popart::TensorId, popart::IArray &> anchors = {{act5, act5_wrapper}};

  WeightsIO weightsRead;
  std::vector<float> w0_readback(weightInfo.nelms(), -99.0f);
  std::vector<float> w1_readback(weightInfo.nelms(), -99.0f);
  std::vector<float> w2_readback(weightInfo.nelms(), -99.0f);
  std::vector<float> w3_readback(weightInfo.nelms(), -99.0f);
  std::vector<float> w4_readback(weightInfo.nelms(), -99.0f);
  std::vector<float> w5_readback(weightInfo.nelms(), -99.0f);
  weightsRead.insert(w0, {w0_readback.data(), weightInfo});
  weightsRead.insert(w1, {w1_readback.data(), weightInfo});
  weightsRead.insert(w2, {w2_readback.data(), weightInfo});
  weightsRead.insert(w3, {w3_readback.data(), weightInfo});
  weightsRead.insert(w4, {w4_readback.data(), weightInfo});
  weightsRead.insert(w5, {w5_readback.data(), weightInfo});

  // write initial weights to host
  session->weightsFromHost();

  float sampleNumVal = 100.0f;
  for (int i = 0; i < 3; ++i) {
    std::cout << "Iteration (call to run(...)) # " << i << std::endl;

    // get new samples
    for (int i = 0; i < samplesPerStep; ++i) {
      for (int j = 0; j < sampleElms; ++j) {
        auto stepIndex       = i * sampleElms + j;
        v_input_0[stepIndex] = fdis(eng) > 0.5 ? -sampleNumVal : +sampleNumVal;
        v_sample_sum_0[j] += v_input_0[stepIndex];
      }
    }
    popart::NDArrayWrapper<float> input1_wrapper(v_input_0.data(),
                                                 stepDataInfo);
    std::map<popart::TensorId, popart::IArray &> inputs = {
        {input1, input1_wrapper}};
    popart::StepIO stepio(inputs, anchors);

    // process the 400 samples (100 batches), streaming back anchors each batch
    session->run(stepio);
  }

  // read final weights back
  session->weightsToHost();
  session->readWeights(weightsRead);

  std::vector<std::vector<float> *> ws = {&w0_readback,
                                          &w1_readback,
                                          &w2_readback,
                                          &w3_readback,
                                          &w4_readback,
                                          &w5_readback};

  // get sum of absolute differences between computed and expected
  float sumAbsDiff = 0.0;
  for (auto wv_p : ws) {
    for (int i = 0; i < w0_readback.size(); ++i) {
      sumAbsDiff += std::abs((*wv_p)[i] + v_sample_sum_0[i] * learnRate *
                                              lambda / sampleNumVal);
    }
  }
  BOOST_CHECK(sumAbsDiff < 1e-5);
}


text:
/*
 * Copyright (c) 2021 The Foundation for Research on Information Technologies in Society (IT'IS).
 *
 * This file is part of iSEG
 * (see https://github.com/ITISFoundation/osparc-iseg).
 *
 * This software is released under the MIT License.
 *  https://opensource.org/licenses/MIT
 */
#include "Precompiled.h"

#include "XdmfImageReader.h"

#include "Data/ScopedTimer.h"
#include "Data/Transform.h"

#include "Core/ColorLookupTable.h"
#include "Core/HDF5Reader.h"

#include <boost/algorithm/string/replace.hpp>

#include <QDir>
#include <QDomDocument>
#include <QDomElement>
#include <QDomText>
#include <QFile>
#include <QFileInfo>
#include <QString>
#include <QTextStream>

#include <vtkSmartPointer.h>

#include <cassert>
#include <stdexcept>
#include <vector>

namespace iseg {

namespace {
int _Read(HDF5Reader& reader, size_t NumberOfSlices, size_t Width, size_t Height, bool ReadContiguousMemory, const std::string& source_dname, float** ImageSlices, const std::string& target_dname, float** WorkSlices, const std::string& tissue_dname, tissues_size_t** TissueSlices)
{
	if (!reader.Exists(source_dname))
	{
		ISEG_WARNING_MSG("no Source array, will initialize to 0...");
		for (int k = 0; k < NumberOfSlices; k++)
		{
			size_t pos = 0;
			for (int j = 0; j < Height; j++)
			{
				for (int i = 0; i < Width; i++, pos++)
				{
					ImageSlices[k][pos] = 0.0f;
				}
			}
		}
	}
	if (!reader.Exists(target_dname))
	{
		ISEG_WARNING_MSG("no Target array, will initialize to 0...");
		for (int k = 0; k < NumberOfSlices; k++)
		{
			size_t pos = 0;
			for (int j = 0; j < Height; j++)
			{
				for (int i = 0; i < Width; i++, pos++)
				{
					WorkSlices[k][pos] = 0.0f;
				}
			}
		}
	}
	if (!reader.Exists(tissue_dname))
	{
		ISEG_WARNING_MSG("no Tissue array, will initialize to 0...");
		for (int k = 0; k < NumberOfSlices; k++)
		{
			size_t pos = 0;
			for (int j = 0; j < Height; j++)
			{
				for (int i = 0; i < Width; i++, pos++)
				{
					TissueSlices[k][pos] = 0;
				}
			}
		}
	}

	if (ReadContiguousMemory)
	{
		// allocate
		size_t const n = NumberOfSlices * Width * Height;
		std::vector<float> buffer_float;
		try
		{
			ISEG_DEBUG("N = " << n << ", bufferFloat.max_size() = " << buffer_float.max_size());
			buffer_float.resize(n);
		}
		catch (std::length_error& le)
		{
			std::cerr << "bufferFloat length error: " << le.what();
			return 0;
		}

		// Source
		if (reader.Exists(source_dname))
		{
			ScopedTimer timer("Read Source");
			if (!reader.Read(&buffer_float[0], source_dname))
			{
				ISEG_ERROR_MSG("reading Source dataset...");
			}
			else
			{
				size_t n = 0;
				for (int k = 0; k < NumberOfSlices; k++)
				{
					size_t pos = 0;
					for (int j = 0; j < Height; j++)
					{
						for (int i = 0; i < Width; i++, pos++)
						{
							ImageSlices[k][pos] = buffer_float[n];
							n++;
						}
					}
				}
			}
		}

		// Target
		if (reader.Exists(target_dname))
		{
			ScopedTimer timer("Read Target");
			if (!reader.Read(&buffer_float[0], target_dname))
			{
				ISEG_ERROR_MSG("reading Target dataset...");
			}
			else
			{
				size_t n = 0;
				for (int k = 0; k < NumberOfSlices; k++)
				{
					size_t pos = 0;
					for (int j = 0; j < Height; j++)
					{
						for (int i = 0; i < Width; i++, pos++)
						{
							WorkSlices[k][pos] = buffer_float[n];
							n++;
						}
					}
				}
			}
		}
		buffer_float.clear();

		// Tissue
		if (reader.Exists(tissue_dname))
		{
			ScopedTimer timer("Read Tissue");
			std::string type;
			std::vector<HDF5Reader::size_type> dims;
			if (!reader.GetDatasetInfo(type, dims, tissue_dname))
			{
				ISEG_ERROR_MSG("reading Tissue data type...");
				return 0;
			}

			if (type == "unsigned char")
			{
				std::vector<unsigned char> buffer_u_char;
				try
				{
					// vector throws a length_error if resized above max_size
					// vector<float> bufferFloat(N);
					ISEG_DEBUG("N = " << n << ", bufferUChar.max_size() = " << buffer_u_char.max_size());
					buffer_u_char.resize(n);
				}
				catch (std::length_error& le)
				{
					ISEG_ERROR("bufferUChar length error: " << le.what());
					return 0;
				}
				if (!reader.Read(&buffer_u_char[0], tissue_dname))
				{
					ISEG_ERROR_MSG("reading Tissue dataset...");
					return 0;
				}
				size_t n = 0;
				for (int k = 0; k < NumberOfSlices; k++)
				{
					size_t pos = 0;
					for (int j = 0; j < Height; j++)
					{
						for (int i = 0; i < Width; i++, pos++)
						{
							TissueSlices[k][pos] =
									(tissues_size_t)buffer_u_char[n];
							n++;
						}
					}
				}
				buffer_u_char.clear();
			}
			else if (type == "unsigned short")
			{
				std::vector<unsigned short> buffer_u_short;
				try
				{
					// vector throws a length_error if resized above max_size
					// vector<float> bufferFloat(N);
					ISEG_DEBUG("N = " << n << ", bufferUShort.max_size() = " << buffer_u_short.max_size());
					buffer_u_short.resize(n);
				}
				catch (std::length_error& le)
				{
					ISEG_ERROR("bufferUShort length error: " << le.what());
					return 0;
				}
				if (!reader.Read(&buffer_u_short[0], tissue_dname))
				{
					ISEG_ERROR_MSG("reading Tissue dataset...");
					return 0;
				}
				size_t n = 0;
				for (int k = 0; k < NumberOfSlices; k++)
				{
					size_t pos = 0;
					for (int j = 0; j < Height; j++)
					{
						for (int i = 0; i < Width; i++, pos++)
						{
							TissueSlices[k][pos] =
									(tissues_size_t)buffer_u_short[n];
							n++;
						}
					}
				}
				buffer_u_short.clear();
			}
			else if (type == "unsigned int")
			{
				std::vector<unsigned int> buffer_u_int;
				try
				{
					// vector throws a length_error if resized above max_size
					// vector<float> bufferFloat(N);
					ISEG_DEBUG("N = " << n << ", bufferUInt.max_size() = " << buffer_u_int.max_size());
					buffer_u_int.resize(n);
				}
				catch (std::length_error& le)
				{
					ISEG_ERROR("bufferUInt length error: " << le.what());
					return 0;
				}
				if (!reader.Read(&buffer_u_int[0], tissue_dname))
				{
					ISEG_ERROR_MSG("reading Tissue dataset...");
					return 0;
				}
				size_t n = 0;
				for (int k = 0; k < NumberOfSlices; k++)
				{
					size_t pos = 0;
					for (int j = 0; j < Height; j++)
					{
						for (int i = 0; i < Width; i++, pos++)
						{
							TissueSlices[k][pos] =
									(tissues_size_t)buffer_u_int[n];
							n++;
						}
					}
				}
				buffer_u_int.clear();
			}
			else
			{
				ISEG_ERROR_MSG("Tissue data type not supported...");
				return 0;
			}
		}
	}
	else
	{
		size_t const slice_size = Width * Height;

		if (reader.Exists(source_dname))
		{
			ScopedTimer timer("Read Source");
			bool ok = true;
			size_t offset = 0;
			for (int k = 0; k < NumberOfSlices; k++, offset += slice_size)
			{
				ok = ok && reader.Read(ImageSlices[k], offset, slice_size, source_dname);
			}
		}
		if (reader.Exists(target_dname))
		{
			ScopedTimer timer("Read Target");
			bool ok = true;
			size_t offset = 0;
			for (int k = 0; k < NumberOfSlices; k++, offset += slice_size)
			{
				ok = ok && reader.Read(WorkSlices[k], offset, slice_size, target_dname);
			}
		}
		if (reader.Exists(tissue_dname))
		{
			ScopedTimer timer("Read Tissue");
			bool ok = true;
			size_t offset = 0;
			for (int k = 0; k < NumberOfSlices; k++, offset += slice_size)
			{
				ok = ok && reader.Read(TissueSlices[k], offset, slice_size, tissue_dname);
			}
		}
	}

	return 1;
}
} // namespace

XdmfImageReader::XdmfImageReader()
{
	this->m_NumberOfSlices = 0;
	this->m_Width = 0;
	this->m_Height = 0;
	this->m_Compression = 1;
	this->m_ImageTransform.SetIdentity();
	this->m_ImageSlices = nullptr;
	this->m_WorkSlices = nullptr;
	this->m_TissueSlices = nullptr;
	this->m_FileName = nullptr;
	this->m_ReadContiguousMemory = false;
}

XdmfImageReader::~XdmfImageReader() { delete[] this->m_FileName; }

int XdmfImageReader::ParseXML()
{
	QFile file(m_FileName);
	if (!file.open(QIODevice::ReadOnly))
	{
		ISEG_ERROR("cannot open " << m_FileName);
		return 0;
	}

	QDomDocument input_document("Slices");
	QByteArray input_content = file.readAll();
	file.close();

	QString msg;
	if (!input_document.setContent(input_content, false, &msg)) // This function is not reentrant.
	{
		ISEG_ERROR("assigning content of " << m_FileName << ": " << msg.toStdString());
	}

	QDomElement root = input_document.documentElement();
	if (root.tagName() != "Xdmf")
	{
		ISEG_ERROR("invalid root element: " << root.tagName().toStdString());
	}

	QDomElement domain = root.firstChildElement(QString("Domain"));
	QDomNode n = domain.firstChild();
	while (!n.isNull())
	{
		QDomElement e = n.toElement();
		if (e.tagName() == "Grid")
		{
			QString type = e.attribute(QString("Type"));
			if (type.isEmpty())
				type = e.attribute(QString("GridType"));
			if (type != QString("Uniform"))
			{
				ISEG_ERROR("unsupported grid type: " << type.toStdString());
				return 0;
			}

			// Geometry
			QDomElement geometry = e.firstChildElement(QString("Geometry"));
			if (geometry.attribute(QString("Type")) !=
							QString("ORIGIN_DXDYDZ") &&
					geometry.attribute(QString("GeometryType")) !=
							QString("ORIGIN_DXDYDZ"))
			{
				ISEG_ERROR_MSG("unsupported geometry type...");
				return 0;
			}
			QDomNodeList list2 = geometry.elementsByTagName("DataItem");
			for (int i2 = 0; i2 < list2.size(); ++i2)
			{
				const QDomElement e2 = list2.item(i2).toElement();
				QString name = e2.attribute(QString("Name"));
				if (name == QString("Origin"))
				{
					QString text = e2.text().trimmed();
					QStringList text_list = text.split(QString(" "));
					if (text_list.size() != 3)
					{
						ISEG_ERROR_MSG("invalid origin...");
						return 0;
					}
					float offset[3];
					for (int i3 = 0; i3 < text_list.size(); ++i3)
						offset[2 - i3] = text_list[i3].toFloat();
					m_ImageTransform.SetOffset(offset);
				}
				else if (name == QString("Spacing"))
				{
					QString text = e2.text().trimmed();
					QStringList text_list = text.split(QString(" "));
					if (text_list.size() != 3)
					{
						ISEG_ERROR_MSG("invalid spacing...");
						return 0;
					}
					for (int i3 = 0; i3 < text_list.size(); ++i3)
						m_PixelSize[2 - i3] = text_list[i3].toFloat();
				}
				else
				{
					ISEG_ERROR_MSG("expecting origin and spacing...");
					return 0;
				}
			}

			// Topology
			QDomElement topology = e.firstChildElement(QString("Topology"));
			if (topology.attribute(QString("Type")) !=
							QString("3DCORECTMesh") &&
					topology.attribute(QString("TopologyType")) !=
							QString("3DCORECTMesh"))
			{
				ISEG_ERROR_MSG("unsupported topology type...");
				return 0;
			}
			QStringList text_list =
					topology.attribute(QString("Dimensions")).split(QString(" "));
			if (text_list.size() != 3)
			{
				ISEG_ERROR_MSG("invalid dimensions...");
				return 0;
			}
			this->m_Width = text_list[2].toInt();
			this->m_Height = text_list[1].toInt();
			this->m_NumberOfSlices = text_list[0].toInt();

			// Attributes
			QDomNodeList attribute_list = e.elementsByTagName("Attribute");
			for (int i2 = 0; i2 < attribute_list.size(); ++i2)
			{
				const QDomElement scalar = attribute_list.item(i2).toElement();
				QString name = scalar.attribute(QString("Name"));
				this->m_ArrayNames.append(name);
				QDomElement dataitem =
						scalar.firstChildElement(QString("DataItem"));
				QString text = dataitem.text().trimmed();
				int idx = text.indexOf(":");
				QString datasetname = text.remove(0, idx + 1);
				this->m_MapArrayNames[name] = datasetname;
			}

		} // if(e.tagName()=="Grid")

		n = n.nextSibling();
	}

	return 1;
}

// BL TODO read rotation & origin from HDF5 part, if available
int XdmfImageReader::Read()
{
	ISEG_INFO("Reading " << this->m_FileName << ": " << m_Width << " x " << m_Height << " x " << m_NumberOfSlices);

	QString q_file_name(this->m_FileName);
	QFileInfo file_info(q_file_name);
	QString basename = file_info.completeBaseName();

	// save working directory
	QDir oldcwd = QDir::current();

	// enter the xmf file folder so relative names for hdf5 files work
	QDir::setCurrent(file_info.absolutePath());

	HDF5Reader reader;
	const QString fname = basename + ".h5";
	if (!reader.Open(fname.toStdString()))
	{
		ISEG_ERROR("opening " << fname.toStdString());
		return 0;
	}

	int r = _Read(reader, m_NumberOfSlices, m_Width, m_Height, m_ReadContiguousMemory, this->m_MapArrayNames["Source"].toAscii().data(), m_ImageSlices, this->m_MapArrayNames["Target"].toAscii().data(), m_WorkSlices, this->m_MapArrayNames["Tissue"].toAscii().data(), m_TissueSlices);

	reader.Close();

	// restore working directory
	QDir::setCurrent(oldcwd.absolutePath());

	return r;
}

std::shared_ptr<ColorLookupTable> XdmfImageReader::ReadColorLookup() const
{
	std::string fname(this->m_FileName);
	boost::replace_last(fname, ".xmf", ".h5");

	HDFImageReader reader;
	reader.SetFileName(fname.c_str());

	return reader.ReadColorLookup();
}

HDFImageReader::HDFImageReader()
{
	this->m_NumberOfSlices = 0;
	this->m_Width = 0;
	this->m_Height = 0;
	this->m_Compression = 1;
	this->m_ImageSlices = nullptr;
	this->m_WorkSlices = nullptr;
	this->m_TissueSlices = nullptr;
	this->m_FileName = nullptr;
	this->m_ReadContiguousMemory = false;
}

HDFImageReader::~HDFImageReader() { delete[] this->m_FileName; }

int HDFImageReader::ParseHDF()
{
	std::cerr << "parsing file " << this->m_FileName << endl;
	HDF5Reader reader;
	if (!reader.Open(m_FileName))
	{
		ISEG_ERROR("opening " << m_FileName);
		return 0;
	}

	bool ok = true;
	std::vector<int> dimensions;
	dimensions.resize(3);
	ok = ok && (reader.Read(dimensions, "/dimensions") != 0) &&
			 (dimensions.size() == 3);
	this->m_Width = dimensions[0];
	this->m_Height = dimensions[1];
	this->m_NumberOfSlices = dimensions[2];

	std::vector<float> pixelsize(3);
	ok = ok && (reader.Read(pixelsize, "/pixelsize") != 0) &&
			 (pixelsize.size() == 3);
	for (unsigned short i = 0; i < 3; i++)
		m_PixelSize[i] = pixelsize[i];

	std::vector<float> offset(3);
	ok = ok && (reader.Read(offset, "/offset") != 0) && (offset.size() == 3);

	if (reader.Exists("/rotation"))
	{
		std::vector<float> rotation(9);
		ok = ok && (reader.Read(rotation, "/rotation") != 0) &&
				 (rotation.size() == 9);

		m_ImageTransform.SetIdentity();
		m_ImageTransform.SetOffset(offset.data());

		for (int k = 0; k < 3; ++k)
		{
			m_ImageTransform[k][0] = rotation[k * 3 + 0];
			m_ImageTransform[k][1] = rotation[k * 3 + 1];
			m_ImageTransform[k][2] = rotation[k * 3 + 2];
		}
	}
	else
	{
		std::vector<float> dc(6);
		ok = ok && (reader.Read(dc, "/dc") != 0) && (dc.size() == 6);

		m_ImageTransform.SetTransform(offset.data(), dc.data());
	}

	//Transform

	if (reader.Exists("/Source"))
	{
		m_MapArrayNames["Source"] = QString("/Source");
		m_ArrayNames.append(QString("Source"));
	}
	if (reader.Exists("/Target"))
	{
		m_MapArrayNames["Target"] = QString("/Target");
		m_ArrayNames.append(QString("Target"));
	}
	if (reader.Exists("/Tissue"))
	{
		m_MapArrayNames["Tissue"] = QString("/Tissue");
		m_ArrayNames.append(QString("Tissue"));
	}

	reader.Close();

	assert(ok);
	return ok ? 1 : 0;
}

int HDFImageReader::Read()
{
	ISEG_INFO("Reading " << this->m_FileName << ": " << m_Width << " x " << m_Height << " x " << m_NumberOfSlices);

	//	vector<int> dims;

	QString q_file_name(this->m_FileName);
	QFileInfo file_info(q_file_name);
	QString basename = file_info.completeBaseName();
	QString suffix = file_info.suffix();

	// save working directory
	QDir oldcwd = QDir::current();

	// enter the xmf file folder so relative names for hdf5 files work
	QDir::setCurrent(file_info.absolutePath());

	HDF5Reader reader;
	const QString fname = basename + "." + suffix;
	if (!reader.Open(fname.toAscii().data()))
	{
		ISEG_ERROR("opening " << fname.toStdString());
		return 0;
	}

	_Read(reader, m_NumberOfSlices, m_Width, m_Height, m_ReadContiguousMemory, this->m_MapArrayNames["Source"].toAscii().data(), m_ImageSlices, this->m_MapArrayNames["Target"].toAscii().data(), m_WorkSlices, this->m_MapArrayNames["Tissue"].toAscii().data(), m_TissueSlices);

	reader.Close();

	// restore working directory
	QDir::setCurrent(oldcwd.absolutePath());

	return 1;
}

std::shared_ptr<ColorLookupTable> HDFImageReader::ReadColorLookup() const
{
	ScopedTimer timer("ReadColorLookup");

	std::shared_ptr<ColorLookupTable> color_lookup_table;

	QString q_file_name(this->m_FileName);
	QFileInfo file_info(q_file_name);
	QString basename = file_info.completeBaseName();
	QString suffix = file_info.suffix();

	// save working directory
	QDir oldcwd = QDir::current();
	QDir::setCurrent(file_info.absolutePath());

	HDF5Reader reader;
	const QString fname = basename + "." + suffix;
	if (!reader.Open(fname.toStdString()))
	{
		ISEG_ERROR("opening " << fname.toStdString());
		return nullptr;
	}

	bool ok = true;
	if (reader.Exists("/Lut"))
	{
		int version = -1, num_colors = 0;
		ok = ok && (reader.Read(&version, "/Lut/version") != 0);
		ok = ok && (reader.Read(&num_colors, "/Lut/size") != 0);

		if (version == 1)
		{
			std::vector<std::string> colors = reader.GetGroupInfo("/Lut");
			if (!colors.empty())
			{
				int index;
				float float_rgb[3];
				unsigned char rgb[3];

				ok = ok && (colors.size() == num_colors + 2);
				if (!ok)
				{
					ISEG_ERROR_MSG("could not load color lookup table");
					return nullptr;
				}
				if (version > 1)
				{
					ISEG_ERROR_MSG("Error: could not load color lookup table. The file "
												 "format is newer than this iSEG.");
					return nullptr;
				}

				color_lookup_table = std::make_shared<ColorLookupTable>();
				color_lookup_table->SetNumberOfColors(num_colors);

				for (const auto& name : colors)
				{
					if (ok && name.find("color") != std::string::npos)
					{
						std::string const folder_name = "/Lut/" + name;

						ok = ok &&
								 (reader.Read(&index, folder_name + "/index") != 0);
						ok = ok &&
								 (reader.Read(float_rgb, folder_name + "/rgb") != 0);

						rgb[0] = static_cast<unsigned char>(float_rgb[0] * 255.0);
						rgb[1] = static_cast<unsigned char>(float_rgb[1] * 255.0);
						rgb[2] = static_cast<unsigned char>(float_rgb[2] * 255.0);
						color_lookup_table->SetColor(index, rgb);
					}
				}
			}
		}
		else if (version == 2)
		{
			std::vector<float> colors(3 * num_colors);
			ok = ok && (reader.Read(colors.data(), "/Lut/colors") != 0);
			if (ok)
			{
				color_lookup_table = std::make_shared<ColorLookupTable>();
				color_lookup_table->SetNumberOfColors(num_colors);

				unsigned char rgb[3];
				for (int i = 0; i < num_colors; ++i)
				{
					rgb[0] = static_cast<unsigned char>(colors[i * 3 + 0] * 255.0);
					rgb[1] = static_cast<unsigned char>(colors[i * 3 + 1] * 255.0);
					rgb[2] = static_cast<unsigned char>(colors[i * 3 + 2] * 255.0);
					color_lookup_table->SetColor(static_cast<size_t>(i), rgb);
				}
			}
		}
		else if (version == 3)
		{
			std::vector<unsigned char> colors(3 * num_colors);
			ok = ok && (reader.Read(colors.data(), "/Lut/colors") != 0);
			if (ok)
			{
				color_lookup_table = std::make_shared<ColorLookupTable>();
				color_lookup_table->SetNumberOfColors(num_colors);

				for (int i = 0; i < num_colors; ++i)
				{
					color_lookup_table->SetColor(static_cast<size_t>(i), &colors[i * 3]);
				}
			}
		}
		else
		{
			ISEG_ERROR_MSG("color lookup table was written with a newer version of iSEG");
		}
	}

	QDir::setCurrent(oldcwd.absolutePath());

	return color_lookup_table;
}

} // namespace iseg

max_stars_repo_name: ITISFoundation/osparc-iseg

max_issues_repo_name: dyollb/osparc-iseg

max_forks_repo_name: dyollb/osparc-iseg