text:
//
// Build a BFieldManager.
//
//

// Includes from C++
#include <fstream>
#include <iomanip>
#include <iostream>
#include <set>

// Includes from C ( needed for block IO ).
#include <errno.h>
#include <fcntl.h>
#include <string.h>
#include <sys/stat.h>
#include <sys/types.h>
#include <unistd.h>

// Includes from boost
#include <boost/iostreams/filter/bzip2.hpp>
#include <boost/iostreams/filter/gzip.hpp>
#include <boost/regex.hpp>
#include "boost/iostreams/filtering_stream.hpp"

// Framework includes
#include "canvas/Utilities/Exception.h"
#include "cetlib_except/exception.h"
#include "messagefacility/MessageLogger/MessageLogger.h"

// Includes from Mu2e
#include "Offline/BFieldGeom/inc/BFInterpolationStyle.hh"
#include "Offline/BFieldGeom/inc/BFieldConfig.hh"
#include "Offline/BFieldGeom/inc/BFieldManager.hh"
#include "Offline/BFieldGeom/inc/DiskRecord.hh"
#include "Offline/GeneralUtilities/inc/MinMax.hh"
#include "Offline/GeometryService/inc/BFieldManagerMaker.hh"

// CLHEP includes
#include "CLHEP/Units/SystemOfUnits.h"

// CSV reader
#include "Offline/GeneralUtilities/inc/csv.hh"

using namespace std;

namespace mu2e {

    namespace {

        typedef std::vector<std::shared_ptr<BFMap>> MapContainerType;
        // To control printouts from helper functions in this file
        int bfieldVerbosityLevel = 0;

        // Passing a string (not a const string&) potentially more efficient b/c of the copying
        std::string basename(std::string file) {
            std::string::size_type i = file.rfind('.');
            if (i != std::string::npos) {
                file.erase(i);
            }

            i = file.rfind('/');
            if (i != std::string::npos) {
                file.erase(0, i + 1);
            }
            return file;
        }
    }  // namespace

    //
    // Utility function which adds filter to decompress file if neccessary
    //
    void decompressFile(const string& filename, boost::iostreams::filtering_istream& in) {
        boost::regex re_gz("^.*\\.gz\\s*$");
        boost::regex re_bz("^.*\\.bz2\\s*$");
        if (boost::regex_match(filename.c_str(), re_gz)) {
            in.push(boost::iostreams::gzip_decompressor());
        } else if (boost::regex_match(filename.c_str(), re_bz)) {
            in.push(boost::iostreams::bzip2_decompressor());
        }
    }

    BFieldManagerMaker::BFieldManagerMaker(const BFieldConfig& config)
        : _resolveFullPath(), _bfmgr(new BFieldManager()) {
        bfieldVerbosityLevel = config.verbosityLevel();

        // break potential mapTypeList into two vectors... kind of ugly right now.
        // Maybe config should just have two mapTypeLists for inner and outer.
        // Eventually, each map file should just be paired with a mapType.
        if (!config.mapTypeList().empty()) {
            _innerTypes = std::vector<BFMapType>(
                config.mapTypeList().begin(),
                config.mapTypeList().begin() + config.innerMapFiles().size());
            _outerTypes =
                std::vector<BFMapType>(config.mapTypeList().begin() + config.innerMapFiles().size(),
                                       config.mapTypeList().end());
        }

        if (config.mapType() == BFMapType::GMC) {
            // Add the field maps.
            for (unsigned i = 0; i < config.gmcDimensions().size(); ++i) {
                readGMCMap(basename(config.outerMapFiles()[i]),
                           _resolveFullPath(config.outerMapFiles()[i]), config.gmcDimensions()[i],
                           config.scaleFactor(), config.interpolationStyle());
            }

        } else if (config.mapType() == BFMapType::G4BL) {
            loadG4BL(&_bfmgr->innerMaps_, config.innerMapFiles(), config.scaleFactor(),
                     config.interpolationStyle());
            loadG4BL(&_bfmgr->outerMaps_, config.outerMapFiles(), config.scaleFactor(),
                     config.interpolationStyle());

        } else if (config.mapType() == BFMapType::PARAM) {
            loadParam(&_bfmgr->innerMaps_, config.innerMapFiles(), _innerTypes,
                      config.interpolationStyle(), config.scaleFactor());
            loadParam(&_bfmgr->outerMaps_, config.outerMapFiles(), _outerTypes,
                      config.interpolationStyle(), config.scaleFactor());

        } else {
            throw cet::exception("GEOM")
                << "Unknown format of file with magnetic field maps: " << config.mapType() << "\n";
        }

        _bfmgr->cm_.setMaps((const MapContainerType&)_bfmgr->innerMaps_,
                            (const MapContainerType&)_bfmgr->outerMaps_);

        // The field manager is fully initialized.
        // Some extra stuff that is convenient to do here:
        if (config.flipBFieldMaps()) {
            for (BFieldManager::MapContainerType::iterator i = _bfmgr->getInnerMaps().begin();
                 i != _bfmgr->getInnerMaps().end(); ++i) {
                flipMap(dynamic_cast<BFGridMap&>(**i));
            }

            for (BFieldManager::MapContainerType::iterator i = _bfmgr->getOuterMaps().begin();
                 i != _bfmgr->getOuterMaps().end(); ++i) {
                flipMap(dynamic_cast<BFGridMap&>(**i));
            }
        }

        if (config.writeBinaries()) {
            for (BFieldManager::MapContainerType::const_iterator i = _bfmgr->getInnerMaps().begin();
                 i != _bfmgr->getInnerMaps().end(); ++i) {
                writeG4BLBinary(dynamic_cast<const BFGridMap&>(**i), (*i)->getKey() + ".bin");
            }

            for (BFieldManager::MapContainerType::const_iterator i = _bfmgr->getOuterMaps().begin();
                 i != _bfmgr->getOuterMaps().end(); ++i) {
                writeG4BLBinary(dynamic_cast<const BFGridMap&>(**i), (*i)->getKey() + ".bin");
            }
        }

        // For debug purposes: print the field in the target region
        if (bfieldVerbosityLevel > 0) {
            CLHEP::Hep3Vector b = _bfmgr->getBField(CLHEP::Hep3Vector(3900.0, 0.0, -6550.0));
            cout << "B-field at the proton target: (" << b.x() << "," << b.y() << "," << b.z()
                 << ")" << endl;
        }
    }

    // Anonymous namespace to hold some helper functions
    namespace {

        // Helper function for parseG4BLHeader.
        void fillGrid(boost::cmatch const& matches,
                      std::vector<double>& X0,
                      std::vector<int>& dim,
                      std::vector<double>& dX) {
            for (int i = 1; i <= 3; i++) {
                double value;
                istringstream sin(string(matches[i].first, matches[i].second));
                sin >> value;
                X0.push_back(value);
            }

            for (int i = 4; i <= 6; i++) {
                int value;
                istringstream sin(string(matches[i].first, matches[i].second));
                sin >> value;
                dim.push_back(value);
            }

            for (int i = 7; i <= 9; i++) {
                double value;
                istringstream sin(string(matches[i].first, matches[i].second));
                sin >> value;
                dX.push_back(value);
            }

        }  // end fillGrid

        // Another helper function for parseG4BLHeader.
        void fillOffset(boost::cmatch const& matches, CLHEP::Hep3Vector& offset) {
            for (int i = 1; i <= 3; i++) {
                double value;
                istringstream sin(string(matches[i].first, matches[i].second));
                sin >> value;
                offset(i - 1) = value;
            }

        }  // end fillOffset

        // Helper function for loadG4BL
        // Parse the header section of a G4BL format .txt or .header file.
        // First argument is an input argument, the others are output arguments.
        void parseG4BLHeader(std::string const& path,
                             std::vector<double>& X0,
                             std::vector<int>& dim,
                             std::vector<double>& dX,
                             CLHEP::Hep3Vector& offset,
                             bool& extendYFound) {
            // The offset parameter is not present in files earlier than Mau7.
            // The value set here is the correct value for all G4BL format files earlier than
            // Mau7. For Mau7 and later, the file will contain a structured comment that
            // contains the correct value.  If the structured comment is present in the file,
            // let it define value of the offset;  if not, use the default given here.
            static const CLHEP::Hep3Vector offsetDefault(-3904.0, 0.0, 7929.0);

            // Open the input file.
            if (path.empty())
                throw "BFieldManagerMaker::loadG4BL: find_file failure!";  // TODO: improve
                                                                           // exception
            ifstream fin(path.c_str());
            if (!fin.is_open()) {
                throw cet::exception("GEOM")
                    << "Could not open file containing the magnetic field data. "
                    << "Filename: " << path << "\n";
            }
            boost::iostreams::filtering_istream in;
            decompressFile(path, in);
            in.push(fin);

            // Regex to parse the string that contains the grid description
            boost::regex reGrid(
                "^\\s*grid"
                "\\s+X0=([eE\\d\\-\\+\\.]+)\\s+Y0=([eE\\d\\-\\+\\.]+)\\s+Z0=([eE\\d\\-\\+\\.]+)"
                "\\s+nX=([eE\\d\\-\\+\\.]+)\\s+nY=([eE\\d\\-\\+\\.]+)\\s+nZ=([eE\\d\\-\\+\\.]+)"
                "\\s+dX=([eE\\d\\-\\+\\.]+)\\s+dY=([eE\\d\\-\\+\\.]+)\\s+dZ=([eE\\d\\-\\+\\.]+)"
                ".*$");

            // Regex to parse the string that contains the offset.
            boost::regex reOffset(
                "^\\s*#\\s+Origin\\s+shift\\s+for\\s+Mu2e:"
                "\\s+([eE\\d\\-\\+\\.]+)\\s+([eE\\d\\-\\+\\.]+)\\s+([eE\\d\\-\\+\\.]+)"
                ".*$");

            // Regex to parse the string that contains the command for map extension.
            boost::regex reExtendY(
                "^\\s*extendY\\s+flip=By"
                ".*$");

            // Regex to parse the string that contains possible additional command for map
            // extension.
            boost::regex reExtend("^\\s*extend.*$");

            // Search for the two lines of interest and, if present, extract their information.
            string cbuf;
            boost::cmatch matches;
            bool gridFound(false);
            bool offsetFound(false);
            extendYFound = false;
            bool extendFound(false);
            bool extendError(false);
            int nread = 100;  // Safety - don't read more than 100 lines
            while ((!in.eof()) && (--nread > 0)) {
                getline(in, cbuf);
                if (boost::regex_match(cbuf.c_str(), matches, reGrid)) {
                    gridFound = true;
                    fillGrid(matches, X0, dim, dX);
                }
                if (boost::regex_match(cbuf.c_str(), matches, reOffset)) {
                    offsetFound = true;
                    fillOffset(matches, offset);
                }
                if (boost::regex_match(cbuf.c_str(), matches, reExtendY)) {
                    extendYFound = true;
                }
                if (boost::regex_match(cbuf.c_str(), matches, reExtend)) {
                    if (!extendYFound || extendFound)
                        extendError = true;
                    extendFound = true;
                }
            }
            in.pop();
            fin.close();

            // The grid description must be present.
            if (!gridFound) {
                throw cet::exception("GEOM")
                    << "Could not find param string in magnetic firld map. "
                    << "Filename: " << path << ", found " << matches.size() << " items."
                    << "\n";
            }

            // The extend command must be correct and unique
            if (extendError) {
                throw cet::exception("GEOM")
                    << "Do not know how to extend magnetic field map. "
                    << "Filename: " << path << ". Only a command extendY flip=By can be managed."
                    << "\n";
            }

            // The offset description is optional and has a default.
            if (!offsetFound) {
                offset = offsetDefault;
                if (bfieldVerbosityLevel > 1) {
                    std::cout << "BFieldManagerMaker: G4BLHeader offset missing, using default = "
                              << offset << std::endl;
                }
            } else {
                if (bfieldVerbosityLevel > 1) {
                    std::cout << "BFieldManagerMaker: using offset from G4BLHeader = " << offset
                              << std::endl;
                }
            }

            // Adjust the lower bounds of the grid.
            X0[0] = X0[0] - offset.x();
            X0[1] = X0[1] - offset.y();
            X0[2] = X0[2] - offset.z();

        }  // end parseHeader

        // string split functions for parsing from csv for param fields
        vector<string>& split(const string& s, char delim, vector<string>& elems) {
            stringstream ss(s);
            string item;
            while (getline(ss, item, delim)) {
                elems.push_back(item);
            }
            return elems;
        }

        vector<string> split(const string& s, char delim) {
            vector<string> elems;
            split(s, delim, elems);
            return elems;
        }

    }  // end anonymous namespace

    // Loads a sequence of Parametric files
    void BFieldManagerMaker::loadParam(BFieldManager::MapContainerType* mapContainer,
                                       const BFieldConfig::FileSequenceType& files,
                                       std::vector<BFMapType> mapTypeList,
                                       BFInterpolationStyle interpStyle,
                                       double scaleFactor) {
        // typedef BFieldConfig::FileSequenceType::const_iterator Iter;

        // for (Iter i = files.begin(); i != files.end(); ++i) {
        for (unsigned i = 0; i < files.size(); ++i) {
            if (bfieldVerbosityLevel > 0) {
                cout << "Reading " << files[i] << endl;
            }
            const std::string mapkey = basename(files[i]);
            BFMapType indivMapType(BFMapType::PARAM);
            if (!mapTypeList.empty())
                indivMapType = mapTypeList[i];
            if (indivMapType == BFMapType::PARAM) {
                loadParam(mapContainer, mapkey, _resolveFullPath(files[i]), scaleFactor);
            } else {
                loadG4BL(mapContainer, mapkey, _resolveFullPath(files[i]), scaleFactor,
                         interpStyle);
            }
        }
    }

    // Loads a sequence of G4BL files

    void BFieldManagerMaker::loadG4BL(BFieldManager::MapContainerType* mapContainer,
                                      const BFieldConfig::FileSequenceType& files,
                                      double scaleFactor,
                                      BFInterpolationStyle interpStyle) {
        typedef BFieldConfig::FileSequenceType::const_iterator Iter;

        for (Iter i = files.begin(); i != files.end(); ++i) {
            if (bfieldVerbosityLevel > 0) {
                cout << "Reading " << *i << endl;
            }

            const std::string mapkey = basename(*i);

            loadG4BL(mapContainer, mapkey, _resolveFullPath(*i), scaleFactor, interpStyle);
        }
    }

    // Parse the config file to learn about one magnetic field map.
    // Create an empty map and call the code to load the map from the file.
    void BFieldManagerMaker::loadParam(BFieldManager::MapContainerType* mapContainer,
                                       const std::string& key,
                                       const std::string& resolvedFileName,
                                       double scaleFactor) {
        // Create an empty map.
        auto dsmap = _bfmgr->addBFParamMap(mapContainer, key, -4696, -3096, -800, 800, 3500, 14000,
                                           BFMapType::PARAM, scaleFactor);
        // Fill the map from the disk file.
        readParamFile(resolvedFileName, *dsmap);
    }

    // Parse the config file to learn about one magnetic field map.
    // Create an empty map and call the code to load the map from the file.
    void BFieldManagerMaker::loadG4BL(BFieldManager::MapContainerType* mapContainer,
                                      const std::string& key,
                                      const std::string& resolvedFileName,
                                      double scaleFactor,
                                      BFInterpolationStyle interpStyle) {
        // Extract information from the header.
        vector<double> X0;
        vector<int> dim;
        vector<double> dX;
        CLHEP::Hep3Vector G4BL_offset;
        bool extendYFound;
        parseG4BLHeader(resolvedFileName, X0, dim, dX, G4BL_offset, extendYFound);

        // Create an empty map.
        BFInterpolationStyle meco(BFInterpolationStyle::meco);
        auto dsmap =
            _bfmgr->addBFGridMap(mapContainer, key, dim[0], X0[0], dX[0], dim[1], X0[1], dX[1],
                                 dim[2], X0[2], dX[2], BFMapType::G4BL, scaleFactor, interpStyle);
        dsmap->_flipy = extendYFound;
        // Fill the map from the disk file.
        if (resolvedFileName.find(".header") != string::npos) {
            readG4BLBinary(resolvedFileName, *dsmap);
        } else {
            readG4BLMap(resolvedFileName, *dsmap, G4BL_offset);
        }
    }

    //
    // Read one magnetic field map file in MECO GMC format.
    //
    // This does a 2 pass operation"
    // 1) Pass 1:
    //      Read the input file into a temporary image in memory.
    //      Find the min and max values of the grid points.
    //      A the end of this pass, compute the grid spacing.
    // 2) Pass 2:
    //      Fill the 3D array from the image in memory.
    //

    void BFieldManagerMaker::readGMCMap(const std::string& mapKey,
                                        const std::string& filename,
                                        const std::vector<int>& dim,
                                        double scaleFactor,
                                        BFInterpolationStyle interpStyle) {
        assert(dim.size() == 3);

        // Open the input file.
        int fd = open(filename.c_str(), O_RDONLY);
        if (!fd) {
            throw cet::exception("GEOM")
                << "Could not open file containing the magnetic filed map for: " << mapKey << "\n"
                << "Filename: " << filename << "\n";
        }

        // Compute number of records in the input file.
        const int nrecords = computeArraySize(fd, filename);

        // Image of the file in memory.
        vector<DiskRecord> data(nrecords, DiskRecord());

        // Read file into memory.
        const int nbytes = nrecords * sizeof(DiskRecord);
        ssize_t s = read(fd, &data[0], nbytes);
        if (s != nbytes) {
            if (s == -1) {
                throw cet::exception("GEOM")
                    << "Error reading magnetic field map: " << mapKey << "\n"
                    << "Filename: " << filename << "\n";
            } else {
                throw cet::exception("GEOM")
                    << "Wrong number of bytes read from magnetic field map: " << mapKey << "\n"
                    << "Filename: " << filename << "\n";
            }
        }

        // Tool to find min and max values of grid points.
        MinMax mmX, mmY, mmZ;

        // Collect distinct values of (X,Y,Z) on the grid points.
        set<float> X, Y, Z;

        // Multiply by this factor to convert from kilogauss to tesla.
        double ratio = CLHEP::kilogauss / CLHEP::tesla;

        // For the image in memory:
        //   1) Transform into the correct set of units.
        //   2) Find min/max of each dimension.
        //   3) Collect unique values of (X,Y,Z) of the grid points.
        for (vector<DiskRecord>::iterator i = data.begin(); i != data.end(); ++i) {
            // Modify in place.
            DiskRecord& r = *i;

            // Unit conversion: from (cm, kG) to (mm,T).
            r.x *= CLHEP::cm;
            r.y *= CLHEP::cm;
            r.z *= CLHEP::cm;
            r.bx *= ratio;
            r.by *= ratio;
            r.bz *= ratio;

            // The one check I can do.
            if (r.head != r.tail) {
                throw cet::exception("GEOM")
                    << "Error reading magnetic field map.  "
                    << "Mismatched head and tail byte counts at record: " << data.size() << "\n"
                    << "Could not open file containing the magnetic filed map for: " << mapKey
                    << "\n"
                    << "Filename: " << filename << "\n";
            }

            // Update min/max information.
            mmX.accumulate(r.x);
            mmY.accumulate(r.y);
            mmZ.accumulate(r.z);

            // Populate the set of all unique grid values.
            X.insert(r.x);
            Y.insert(r.y);
            Z.insert(r.z);
        }

        // Expected grid dimentsions.
        const size_t nx = dim[0];
        const size_t ny = dim[1];
        const size_t nz = dim[2];

        // Cross-check that the grid read from the file has the size we expected.
        // This is not really a perfect check since there could be round off error
        // in the computation of the grid points.  But the MECO GMC files were written
        // in a way that this test works.
        if (X.size() != nx || Y.size() != ny || Z.size() != nz) {
            throw cet::exception("GEOM")
                << "Mismatch in expected and observed number of grid points for BField map: "
                << mapKey << "\n"
                << "From file: " << filename << "\n"
                << "Expected/Observed x: " << nx << " " << X.size() << "\n"
                << "Expected/Observed y: " << ny << " " << Y.size() << "\n"
                << "Expected/Observed z: " << nz << " " << Z.size() << "\n";
        }

        // Cross-check that we did not find more data than we have room for.
        if (data.size() > nx * ny * nz - 1) {
            throw cet::exception("GEOM")
                << "Too many values read into the field map: " << mapKey << "\n"
                << "From file: " << filename << "\n"
                << "Expected/Observed size: " << nx * ny * nz << " " << data.size() << "\n";
        }

        // Create an empty map.
        auto bfmap = _bfmgr->addBFGridMap(
            &_bfmgr->outerMaps_, mapKey, nx, mmX.min(), (mmX.max() - mmX.min()) / (nx - 1), ny,
            mmY.min(), (mmY.max() - mmY.min()) / (ny - 1), nz, mmZ.min(),
            (mmZ.max() - mmZ.min()) / (nz - 1), BFMapType::GMC, scaleFactor, interpStyle);

        // Store grid points and field values into 3D arrays
        for (vector<DiskRecord>::const_iterator i = data.begin(), e = data.end(); i != e; ++i) {
            DiskRecord const& r(*i);

            // Find indices corresponding to this grid point.
            // By construction the indices must be in bounds ( we set the limits above ).
            std::size_t ix = bfmap->iX(r.x);
            std::size_t iy = bfmap->iY(r.y);
            std::size_t iz = bfmap->iZ(r.z);

            // Store the information into the 3d arrays.
            bfmap->_field.set(ix, iy, iz, CLHEP::Hep3Vector(r.bx, r.by, r.bz));
            bfmap->_isDefined.set(ix, iy, iz, true);
        }

        return;
    }

    //
    // Read one magnetic field map file in G4BL (TD) format.
    //
    void BFieldManagerMaker::readG4BLMap(const string& filename,
                                         BFGridMap& bfmap,
                                         CLHEP::Hep3Vector G4BL_offset) {
        // Open the input file.
        ifstream fin(filename.c_str());
        if (!fin.is_open())
            throw cet::exception("GEOM") << "Could not open file " << filename << "\n";

        // Skip lines until "data" keyword
        boost::iostreams::filtering_istream in;
        decompressFile(filename, in);
        in.push(fin);
        std::string cbuf;
        boost::regex re("^\\s*data.*$");
        while (!in.eof()) {
            getline(in, cbuf);
            if (boost::regex_match(cbuf, re))
                break;
        }
        if (in.eof())
            throw cet::exception("GEOM") << "Can't find data keyword in " << filename << "\n";

        // Expected grid dimentsions.
        const int nx = bfmap._nx;
        const int ny = bfmap._ny;
        const int nz = bfmap._nz;

        // Calculate expected number of lines to read
        const int nrecord = nx * ny * nz;

        // Read data
        double x[3], b[3];
        int nread = 0;
	getline(in, cbuf);
        while (!in.eof()) {

            istringstream sin(cbuf);
            if ((sin >> x[0] >> x[1] >> x[2] >> b[0] >> b[1] >> b[2]).fail()) {
                break;
		throw cet::exception("BFIELD_FILE_CORRUPT") << "Could not parse  " << filename << "\nat line: "<<cbuf<<"\n" ;
	    }

            // sucessfully read and parsed a line
            ++nread;

            if (nread > nrecord) {
	      throw cet::exception("BFIELD_FILE_TOO_LONG") << "Found too many lines in " << filename << "\n";
	    }

            // Store the information into the 3d arrays.
            CLHEP::Hep3Vector pos = CLHEP::Hep3Vector(x[0], x[1], x[2]) - G4BL_offset;

            // find out the closest grid point
            BFGridMap::GridPoint ipos(bfmap.point2grid(pos));

            // check the range
            if (!bfmap.isValid(ipos)) {
                throw cet::exception("GEOM")
                    << "Error while reading " << filename << "\n"
                    << "line: " << cbuf << "\n"
                    << "the coordinates are out of range: grid point would be " << ipos;
            }

            // Make sure the coordinates we read match a grid point position
            const double tolerance = std::numeric_limits<double>::epsilon();
            CLHEP::Hep3Vector cellfrac(bfmap.cellFraction(pos, ipos));
            if ((std::abs(cellfrac.x()) > tolerance) || (std::abs(cellfrac.y()) > tolerance) ||
                (std::abs(cellfrac.z()) > tolerance)) {
                throw cet::exception("GEOM") << "Error while reading " << filename << "\n"
                                             << "line: " << cbuf << "\n"
                                             << "the coordinates do not match position of any grid "
                                                "point within the tolerance="
                                             << tolerance;
            }

            if (bfmap._isDefined(ipos.ix, ipos.iy, ipos.iz)) {
                throw cet::exception("GEOM") << "Error while reading " << filename << "\n"
                                             << "line: " << cbuf << "\n"
                                             << "Attempt to redefine an already defined point.";
            }

            bfmap._field.set(ipos.ix, ipos.iy, ipos.iz,
			    CLHEP::Hep3Vector(b[0], b[1], b[2]));
            bfmap._isDefined.set(ipos.ix, ipos.iy, ipos.iz, true);

	    // next line
            getline(in, cbuf);
        }

        if (nread != nrecord) {
            throw cet::exception("GEOM")
                << "Error while reading " << filename << "\n"
                << "Read " << nread << " out of expected " << nrecord << " lines.\n"
                << "Last line:\n"
                << cbuf << "\n";
        }

        return;
    }

    void BFieldManagerMaker::readG4BLBinary(const string& headerFilename, BFGridMap& bf) {
        // Form the name of the binary file from the name of the header file.
        string::size_type i = headerFilename.find(".header");
        if (i == string::npos) {
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:readG4BLBinary Expected a file type of .header: "
                << headerFilename << "\n";
        }
        string binFilename = headerFilename.substr(0, i);
        binFilename += ".bin";  // bin file should come from the same dir as the header - do not
                                // search the path

        // Number of points in each big array.
        int nPoints = bf.nx() * bf.ny() * bf.nz();

        // Number of bytes in each big array.
        size_t nbytes = sizeof(CLHEP::Hep3Vector) * nPoints;

        // Open the binary file.
        int fd = open(binFilename.c_str(), O_RDONLY);
        if (fd < 0) {
            int errsave = errno;
            char* errmsg = strerror(errsave);
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:readG4BLBinary Error opening " << binFilename
                << "  errno: " << errsave << " " << errmsg << "\n";
        }

        // Mu2e binary field map files used to store an array of grid point coordinates.
        // This is not done any more.  Here we figure out whether we are dealing with the
        // old or the new format.

        bool grid_coordinates_stored = false;
        struct stat info;
        if (fstat(fd, &info)) {
            int errsave = errno;
            char* errmsg = strerror(errsave);
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:readG4BLBinary: Error doing fstat() on " << binFilename
                << "  errno: " << errsave << " " << errmsg << "\n";
        }
        if (unsigned(info.st_size) == sizeof(unsigned) + 2 * nbytes) {
            grid_coordinates_stored = true;
        } else {  // make sure the new format file has the expected size
            if (unsigned(info.st_size) != sizeof(unsigned) + nbytes) {
                throw cet::exception("GEOM")
                    << "BFieldManagerMaker:readG4BLBinary(): the size = " << info.st_size
                    << " of the file " << binFilename
                    << " does not match any of the expected values (" << sizeof(unsigned) + nbytes
                    << " or " << sizeof(unsigned) + 2 * nbytes << ")";
            }
        }
        if (bfieldVerbosityLevel > 1) {
            std::cout << "grid_coordinates_stored = " << grid_coordinates_stored << " for file "
                      << binFilename << std::endl;
        }

        // A word to receive the endian marker.
        unsigned int marker(0);

        // Address of the first element in the field array
        CLHEP::Hep3Vector* fieldAddr = &bf._field.get(0, 0, 0);

        // Read the endian marker and check that it is OK.
        ssize_t s0 = read(fd, &marker, sizeof(unsigned int));
        if (s0 != sizeof(unsigned int)) {
            int errsave = errno;
            char* errmsg = strerror(errsave);
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:readG4BLBinary Error reading endian marker from "
                << binFilename << "\n"
                << "Status: " << s0 << "  errno: " << errsave << " " << errmsg << "\n";
        }

        static const unsigned int deadbeef(0XDEADBEEF);
        if (marker != deadbeef) {
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:readG4BLBinary endian mismatch" << binFilename
                << "  returned value: " << std::hex << marker << "  expected value: " << deadbeef
                << std::dec << "\n"
                << "Suggestion: change from binary format field maps to the text or gzipped "
                   "text "
                   "format.\n"
                << "\n";
        }

        // Old mu2e binary files had an array of grid pont coordinates stored
        // Need to skip this to get to the the actual B field data
        if (grid_coordinates_stored) {
            off_t res = lseek(fd, nbytes, SEEK_CUR);
            if (res == off_t(-1)) {
                int errsave = errno;
                char* errmsg = strerror(errsave);
                throw cet::exception("GEOM") << "BFieldManagerMaker:readG4BLBinary(): Error doing "
                                                "lseek() on the old format binary "
                                             << binFilename << "\n"
                                             << "  errno: " << errsave << " " << errmsg << "\n";
            }
        }
        // Read the field values at the grid points.
        size_t s2 = read(fd, fieldAddr, nbytes);
        if (s2 != nbytes) {
            int errsave = errno;
            char* errmsg = strerror(errsave);
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:readG4BLBinary Error reading field values from "
                << binFilename << "\n"
                << "Status: " << s2 << "  errno: " << errsave << " " << errmsg << "\n";
        }

        // These maps fill the full box so mark all grid points as valid.
        for (int ix = 0; ix < bf.nx(); ++ix) {
            for (int iy = 0; iy < bf.ny(); ++iy) {
                for (int iz = 0; iz < bf.nz(); ++iz) {
                    bf._isDefined.set(ix, iy, iz, true);
                }
            }
        }

        close(fd);

    }  // end BFieldManagerMaker::readG4BLBinary

    //
    // Read one magnetic field parameter csv.
    //
    void BFieldManagerMaker::readParamFile(const string& filename, BFParamMap& bfmap) {
        // Open the input file.
        ifstream fin(filename.c_str());
        if (!fin.is_open())
            throw cet::exception("GEOM") << "Could not open file " << filename << "\n";

        // Read in csv, tokenize the parameter name, start filling
        // the private class members with the correct values.
        // This assumes the csv is ordered correctly, so I should have some assertions
        // that raise exceptions when ordering is incorrect.  Read-time is not critical,
        // as this is only read from file once, then stored in memory.
        io::CSVReader<2> in(filename);
        in.set_header("param", "val");
        string param;
        double val(-1);
        int ns(-1);
        int ms(-1);
        double Reff(-1);
        vector<vector<double>> As;
        vector<vector<double>> Bs;
        vector<double> Ds;
        // Grab the first 3 vals
        while (in.read_row(param, val)) {
            vector<string> tparams = split(param, '_');
            if (tparams[0].compare("R") == 0) {
                Reff = val;
            } else if (tparams[0].compare("ns") == 0) {
                ns = val;
            } else if (tparams[0].compare("ms") == 0) {
                ms = val;
            }
            if (!(ns == -1 || ms == -1 || Reff == -1))
                break;
        }
        // Ready the 2D arrays
        for (int i = 0; i < ns; ++i) {
            As.push_back(vector<double>());
            Bs.push_back(vector<double>());
        }
        // Fill the 2D params
        while (in.read_row(param, val)) {
            vector<string> tparams = split(param, '_');
            if (tparams[0].compare("A") == 0) {
                As[stoi(tparams[1])].push_back(val);
            } else if (tparams[0].compare("B") == 0) {
                Bs[stoi(tparams[1])].push_back(val);
            } else if (tparams[0].compare("D") == 0) {
                Ds.push_back(val);
            }
        }

        bfmap._ns = ns;
        bfmap._ms = ms;
        bfmap._Reff = Reff;
        bfmap._As = As;
        bfmap._Bs = Bs;
        bfmap._Ds = Ds;

        bfmap.calcConstants();

        cout << Bs[2][69] << endl;
        cout << Reff << ns << ms << endl;

        return;
    }  // namespace mu2e

    void BFieldManagerMaker::writeG4BLBinary(const BFGridMap& bf, const std::string& outputfile) {
        // Number of points in the big array.
        int nPoints = bf.nx() * bf.ny() * bf.nz();

        // Number of bytes in the big array.
        size_t nBytes = sizeof(CLHEP::Hep3Vector) * nPoints;

        cout << "Writing G4BL Magnetic field map in binary format to file: " << outputfile << endl;

        // A marker to catch endian mismatch on readback.
        unsigned int deadbeef(0XDEADBEEF);

        // Address of the first element in the big array.
        CLHEP::Hep3Vector const* fieldAddr = &bf._field.get(0, 0, 0);

        // Open the output file.
        mode_t mode = S_IRUSR | S_IWUSR | S_IRGRP | S_IROTH;
        int flags = O_CREAT | O_WRONLY | O_TRUNC | O_EXCL;
        int fd = open(outputfile.c_str(), flags, mode);
        int errsave = errno;

        // Check for errors.
        if (fd < 0) {
            if (errsave == EEXIST) {
                throw cet::exception("GEOM") << "BFieldManagerMaker:writeG4BLBinary Error opening "
                                             << outputfile << "  File already exists.\n";
            }
            char* errmsg = strerror(errsave);
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:writeG4BLBinary Error opening " << outputfile
                << "  errno: " << errsave << " " << errmsg << "\n";
        }

        // Write the endian marker and the big arrays.
        ssize_t s0 = write(fd, &deadbeef, sizeof(unsigned int));
        if (s0 == -1) {
            errsave = errno;
            char* errmsg = strerror(errsave);
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:writeG4BLBinary Error writing endian marker to "
                << outputfile << "  errno: " << errsave << " " << errmsg << "\n";
        }
        ssize_t s2 = write(fd, fieldAddr, nBytes);
        if (s2 == -1) {
            errsave = errno;
            char* errmsg = strerror(errsave);
            throw cet::exception("GEOM")
                << "BFieldManagerMaker:writeG4BLBinary Error writing field values to " << outputfile
                << "  errno: " << errsave << " " << errmsg << "\n";
        }

        close(fd);

        cout << "Writing complete for file: " << outputfile << endl;

    }  // end BFieldManagerMaker::writeG4BLBinary

    // Compute the size of the array needed to hold the raw data of the field map.
    int BFieldManagerMaker::computeArraySize(int fd, const string& filename) {
        // Get the file size, in bytes, ( info.st_size ).
        struct stat info;
        if (fstat(fd, &info)) {
            int errsave = errno;
            char* errmsg = strerror(errsave);
            throw cet::exception("GEOM")
                << "BFieldManagerMaker::computeArraySize(): Error doing fstat() on " << filename
                << "  errno: " << errsave << " " << errmsg << "\n";
        }

        // Check that an integral number of records fits in the file.
        int remainder = info.st_size % sizeof(DiskRecord);
        if (remainder != 0) {
            throw cet::exception("GEOM")
                << "Field map file does not hold an integral number of records: \n"
                << "Filename:  " << filename << "\n"
                << "Size:      " << info.st_size << "\n"
                << "Remainder: " << remainder << "\n";
        }

        // Compute number of records.
        int nrecords = info.st_size / sizeof(DiskRecord);
        return nrecords;
    }

    void BFieldManagerMaker::flipMap(BFGridMap& bf) {
        std::cout << "Flipping B field vector in map " << bf.getKey() << std::endl;
        for (int ix = 0; ix < bf.nx(); ++ix) {
            for (int iy = 0; iy < bf.ny(); ++iy) {
                for (int iz = 0; iz < bf.nz(); ++iz) {
                    bf._field.set(ix, iy, iz, -bf._field.get(ix, iy, iz));
                }
            }
        }
    }
}  // end namespace mu2e

ext: cc

max_stars_repo_name: resnegfk/Offline

max_issues_repo_name: resnegfk/Offline

max_forks_repo_name: resnegfk/Offline


text:
/*
 * Software License Agreement (BSD-3 License)
 *
 * Copyright (c) 2018 Daniel Koch.
 * All rights reserved.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions are met:
 *
 * * Redistributions of source code must retain the above copyright notice, this
 *   list of conditions and the following disclaimer.
 *
 * * Redistributions in binary form must reproduce the above copyright notice,
 *   this list of conditions and the following disclaimer in the documentation
 *   and/or other materials provided with the distribution.
 *
 * * Neither the name of the copyright holder nor the names of its
 *   contributors may be used to endorse or promote products derived from
 *   this software without specific prior written permission.
 *
 * 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.
 */

/**
 * @file comm.cpp
 * @author Daniel Koch <danielpkoch@gmail.com>
 */

#include <async_comm/comm.h>

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

namespace async_comm
{

DefaultMessageHandler default_message_handler_;

Comm::Comm(MessageHandler& message_handler) :
  message_handler_(message_handler),
  io_service_(),
  new_data_(false),
  shutdown_requested_(false),
  write_in_progress_(false)
{
}

Comm::~Comm()
{
}

bool Comm::init()
{
  if (!do_init())
    return false;

  callback_thread_ = std::thread(std::bind(&Comm::process_callbacks, this));

  async_read();
  io_thread_ = std::thread(boost::bind(&boost::asio::io_service::run, &this->io_service_));

  return true;
}

void Comm::close()
{
  // send shutdown signal to callback thread
  {
    std::unique_lock<std::mutex> lock(callback_mutex_);
    shutdown_requested_ = true;
  }
  condition_variable_.notify_one();

  io_service_.stop();
  do_close();

  if (io_thread_.joinable())
  {
    io_thread_.join();
  }

  if (callback_thread_.joinable())
  {
    callback_thread_.join();
  }
}

void Comm::send_bytes(const uint8_t *src, size_t len)
{
  mutex_lock lock(write_mutex_);

  for (size_t pos = 0; pos < len; pos += WRITE_BUFFER_SIZE)
  {
    size_t num_bytes = (len - pos) > WRITE_BUFFER_SIZE ? WRITE_BUFFER_SIZE : (len - pos);
    write_queue_.emplace_back(src + pos, num_bytes);
  }

  async_write(true);
}

void Comm::register_receive_callback(std::function<void(const uint8_t*, size_t)> fun)
{
  receive_callback_ = fun;
}

void Comm::register_listener(CommListener &listener)
{
  listeners_.push_back(listener);
}

void Comm::async_read()
{
  if (!is_open()) return;

  do_async_read(boost::asio::buffer(read_buffer_, READ_BUFFER_SIZE),
                boost::bind(&Comm::async_read_end,
                            this,
                            boost::asio::placeholders::error,
                            boost::asio::placeholders::bytes_transferred));
}

void Comm::async_read_end(const boost::system::error_code &error, size_t bytes_transferred)
{
  if (error)
  {
    message_handler_.error(error.message());
    close();
    return;
  }

  {
    std::unique_lock<std::mutex> lock(callback_mutex_);
    read_queue_.emplace_back(read_buffer_, bytes_transferred);
    new_data_ = true;
  }
  condition_variable_.notify_one();

  async_read();
}

void Comm::async_write(bool check_write_state)
{
  if (check_write_state && write_in_progress_)
    return;

  mutex_lock lock(write_mutex_);
  if (write_queue_.empty())
    return;

  write_in_progress_ = true;
  WriteBuffer& buffer = write_queue_.front();
  do_async_write(boost::asio::buffer(buffer.dpos(), buffer.nbytes()),
                 boost::bind(&Comm::async_write_end,
                             this,
                             boost::asio::placeholders::error,
                             boost::asio::placeholders::bytes_transferred));
}

void Comm::async_write_end(const boost::system::error_code &error, size_t bytes_transferred)
{
  if (error)
  {
    message_handler_.error(error.message());
    close();
    return;
  }

  mutex_lock lock(write_mutex_);
  if (write_queue_.empty())
  {
    write_in_progress_ = false;
    return;
  }

  WriteBuffer& buffer = write_queue_.front();
  buffer.pos += bytes_transferred;
  if (buffer.nbytes() == 0)
  {
    write_queue_.pop_front();
  }

  if (write_queue_.empty())
    write_in_progress_ = false;
  else
    async_write(false);
}

void Comm::process_callbacks()
{
  std::list<ReadBuffer> local_queue;

  while (true)
  {
    // wait for either new data or a shutdown request
    std::unique_lock<std::mutex> lock(callback_mutex_);
    condition_variable_.wait(lock, [this]{ return new_data_ || shutdown_requested_; });

    // if shutdown requested, end thread execution
    if (shutdown_requested_)
    {
      break;
    }

    // move data to local buffer
    local_queue.splice(local_queue.end(), read_queue_);

    // release mutex to allow continued asynchronous read operations
    new_data_ = false;
    lock.unlock();

    // execute callbacks for all new data
    while (!local_queue.empty())
    {
      ReadBuffer buffer = local_queue.front();
      if (receive_callback_)
      {
        receive_callback_(buffer.data, buffer.len);
      }
      for (std::reference_wrapper<CommListener> listener_ref : listeners_)
      {
        listener_ref.get().receive_callback(buffer.data, buffer.len);
      }
      local_queue.pop_front();
    }
  }
}

} // namespace async_comm

ext: cpp


text:
/*
Copyright (c) 2012 by Procera Networks, Inc. ("PROCERA")

Permission to use, copy, modify, and/or distribute this software for
any purpose with or without fee is hereby granted, provided that the
above copyright notice and this permission notice appear in all
copies.

THE SOFTWARE IS PROVIDED "AS IS" AND PROCERA DISCLAIMS ALL WARRANTIES
WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL PROCERA BE LIABLE FOR
ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT
OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
*/
#include "koi.hpp"

#define KOI_DAEMONIZE 1

#if KOI_DAEMONIZE
#include <sys/types.h>
#include <fcntl.h>
#endif//KOI_DAEMONIZE

#include <stdlib.h>
#include <stdio.h>
#include <string>
#include <vector>
#include <map>
#include <boost/program_options.hpp>
#include "msg.hpp"
#include "file.hpp"

using namespace std;
using namespace boost;
using namespace koi;
using namespace boost::program_options;

extern const char* g_program_build;

namespace {
    void usage(string const& hlp) {
        printf("koi %08x (%s)\n\n"                    \
               "   >(>     ><>               ><'>\n"    \
               "             ><((\">\n"                 \
               "   ><(((*>\n"                           \
               "%s"                                     \
               "Services should be placed in %s.\n"     \
               "\n", koi::version, g_program_build, hlp.c_str(), KOI_CONFIG_SERVICES);
    }

#if KOI_DAEMONIZE
    void daemonize_now(const char* pidfile) {
        pid_t pid = fork();
        if (pid < 0) // fork failed
            exit(1);
        if (pid > 0) { // parent
            if (strlen(pidfile) != 0) {
                file pf(pidfile, "w");
                if (pf.open()) {
                    fprintf(pf._file, "%d\n", (int)pid);
                }
            }
            exit(0);
        }

        // clear file creation mode mask
        umask(0);

        // new SID for child process
        pid_t sid = setsid();
        if (sid < 0)
            exit(1);

        // set workdir so we don't occupy
        // any temporary filesystems
        if ((chdir("/")) < 0)
            exit(1);

        // don't talk to strangers
        freopen("/dev/null", "r", stdin);
        freopen("/dev/null", "w", stdout);
        freopen("/dev/null", "w", stderr);
    }
#endif//KOI_DAEMONIZE

    bool parse_commandline(int argc, char* argv[], vector<string>& configs, string& name) {
            options_description generic("Generic options");
            generic.add_options()
                ("version,v", "Print version string")
                ("help,h", "Produce help message")
                ("no-color,C", "Disable log color output")
#if KOI_DAEMONIZE
                ("daemonize", value<string>(), "Detach from parent process (arg: PID file)")
#endif
                ("log", value<string>()->default_value("console"), "Log targets (console, syslog)")
                ;

            options_description conf("Configuration options");

            string def_file(KOI_CONFIG_FILES);
            typedef vector<string> strings;

            conf.add_options()
                ("file,f", value<strings>()->default_value(configs, def_file), "Configuration file")
                ("name,n", value<string>(), "Set the name of this node")
                ("debug,d", "Debug mode");

            options_description cmdline_options("");
            cmdline_options.add(generic).add(conf);

            positional_options_description p;

            variables_map vm;
            store(command_line_parser(argc, argv).options(cmdline_options).positional(p).run(), vm);
            notify(vm);

            if (vm.count("help")) {
                stringstream ss;
                ss << cmdline_options << endl;
                usage(ss.str());
                return false;
            }

            if (vm.count("version")) {
                printf("koi %08x (%s)\n", koi::version, g_program_build);
                return false;
            }

            if (vm.count("debug")) {
                koi::debug_mode = true;
            }

            if (vm.count("name")) {
                name = vm["name"].as<string>();
            }

#if KOI_DAEMONIZE
            if (vm.count("daemonize")) {
                string pidfile = vm["daemonize"].as<string>();
                daemonize_now(pidfile.c_str());
            }
#endif

            int logmodes = 0;
            if (vm.count("log")) {
                string logs = vm["log"].as<string>();
                if (logs.find("syslog") != string::npos) {
                    logmodes |= logging::LogToSyslog;
                }
                if (logs.find("console") != string::npos) {
                    logmodes |= logging::LogToConsole;
                }
                if (logs.find("file") != string::npos) {
                    logmodes |= logging::LogToFile;
                }
            }
            else {
                logmodes = logging::LogToConsole;
            }

            logging::set_log_mode(logmodes);
            logging::colors = !vm.count("no-color");

            configs = vm["file"].as<vector<string>>();

            return true;
    }
}

int main(int argc, char* argv[]) {
    try {
        vector<string> configs = split(KOI_CONFIG_FILES, ";");
        string name("");

        if (!parse_commandline(argc, argv, configs, name))
            return 0;

        settings cfg;
        if (name != "")
            cfg._name = name;
        if (!cfg.boot(configs))
            return 1;
        return run(cfg, configs);
    }
    catch(const unknown_option& e) {
        fprintf(stderr, "Unknown option: %s. See -h or --help.\n", e.what());
        return 1;
    }
    catch(const too_many_positional_options_error& e) {
        fprintf(stderr, "Error: Illegal positional arguments. See -h or --help.\n");
        return 1;
    }
    catch(const std::exception& e) {
        fprintf(stderr, "Internal error: %s\n", e.what());
        return 1;
    }
}

max_stars_repo_name: yogarajbaskaravel/koi

max_issues_repo_name: yogarajbaskaravel/koi

max_forks_repo_name: yogarajbaskaravel/koi


text:
/**
* @license Apache-2.0
*
* Copyright (c) 2018 The Stdlib Authors.
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
*    http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/

/**
* Generate Boost test fixtures.
*
* ## Notes
*
* -   Run this script from the directory in which fixtures should be written.
*
*/
#include <stdlib.h>
#include <stdio.h>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_int_distribution.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include <boost/math/special_functions/factorials.hpp>

using boost::random::uniform_real_distribution;
using boost::random::uniform_int_distribution;
using boost::random::mt19937;

// Define a new pseudorandom number generator:
mt19937 rng;

/**
* Generates a linearly spaced numeric array of doubles.
*
* ## Notes
*
* -   Assumes that the output array has at least 2 elements.
* -   Output array is guaranteed to include both the start and end values.
*
*
* @param out    output array
* @param len    array length
* @param start  first value
* @param end    last value
*/
void linspace_f64( double *out, const unsigned int len, const double start, const double end ) {
	unsigned int i;
	double incr;

	incr = (end-start) / (len-1);
	for ( i = 0; i < len-1; i++ ) {
		out[ i ] = start + (incr*i);
	}
	out[ i ] = end;
}

/**
* Generates a linearly spaced numeric array of integers.
*
* ## Notes
*
* -   Assumes that the output array has at least 2 elements.
* -   Output array is guaranteed to include both the start and end values.
*
*
* @param out    output array
* @param len    array length
* @param start  first value
* @param end    last value
*/
void linspace_i32( int *out, const unsigned int len, const int start, const int end ) {
	unsigned int i;
	int incr;

	incr = (end-start) / (len-1);
	for ( i = 0; i < len-1; i++ ) {
		out[ i ] = start + (incr*i);
	}
	out[ i ] = end;
}

/**
* Generates an array of pseudorandom doubles drawn from a uniform distribution.
*
* @param out  output array
* @param len  array length
* @param a    lower bound (inclusive)
* @param b    upper bound (exclusive)
*/
void rand_array_f64( double *out, const unsigned int len, const double a, const double b ) {
	unsigned int i;

	// Define a uniform distribution for generating pseudorandom numbers:
	uniform_real_distribution<> randu( a, b );

	for ( i = 0; i < len; i++ ) {
		out[ i ] = randu( rng );
	}
}

/**
* Generates an array of pseudorandom integers drawn from a uniform distribution.
*
* @param out  output array
* @param len  array length
* @param a    lower bound (inclusive)
* @param b    upper bound (exclusive)
*/
void rand_array_i32( int *out, const unsigned int len, const int a, const int b ) {
	unsigned int i;

	// Define a uniform distribution for generating pseudorandom numbers:
	uniform_int_distribution<> randi( a, b );

	for ( i = 0; i < len; i++ ) {
		out[ i ] = randi( rng );
	}
}

/**
* Writes an array of doubles to a file as a series of comma-separated values.
*
* @param f    file to write to
* @param x    array of doubles
* @param len  array length
*/
void write_array_f64( FILE *f, const double *x, const unsigned int len ) {
	unsigned int i;
	for ( i = 0; i < len; i++ ) {
		fprintf( f, "%.17g", x[ i ] );
		if ( i < len-1 ) {
			fprintf( f, "," );
		}
	}
}

/**
* Writes an array of integers to a file as a series of comma-separated values.
*
* @param f    file to write to
* @param x    array of integers
* @param len  array length
*/
void write_array_i32( FILE *f, const int *x, const unsigned int len ) {
	unsigned int i;

	for ( i = 0; i < len; i++ ) {
		fprintf( f, "%d", x[ i ] );
		if ( i < len-1 ) {
			fprintf( f, "," );
		}
	}
}

/**
* Writes a named array of doubles to a file as JSON.
*
* @param f     file to write to
* @param name  array name
* @param x     data
* @param len   array length
*/
void write_named_array_f64( FILE *f, const char *name, const double *x, const unsigned int len ) {
	fprintf( f, "\"%s\":[", name );
	write_array_f64( f, x, len );
	fprintf( f, "]" );
}

/**
* Writes a named array of integers to a file as JSON.
*
* @param f     file to write to
* @param name  array name
* @param x     data
* @param len   array length
*/
void write_named_array_i32( FILE *f, const char *name, const int *x, const unsigned int len ) {
	fprintf( f, "\"%s\":[", name );
	write_array_i32( f, x, len );
	fprintf( f, "]" );
}

/**
* Writes data to a file as JSON.
*
* ## Notes
*
* -   This function SHOULD be tailored to the input data (e.g., input types, output types, number of arguments, etc) and may vary from use case to use case.
*
*
* @param f    file to write to
* @param x    first parameter
* @param n    second parameter
* @param y    results
* @param len  array length
*/
void write_data_as_json( FILE *f, const double *x, const int *n,const double *y, const unsigned int len ) {
	fprintf( f, "{" );
	write_named_array_f64( f, "x", x, len );
	fprintf( f, "," );
	write_named_array_i32( f, "n", n, len );
	fprintf( f, "," );
	write_named_array_f64( f, "expected", y, len );
	fprintf( f, "}" );
}

/**
* Generates test fixtures.
*
* @param x     domain
* @param len   number of values in the domain
* @param name  output filename
*/
void generate( double *x, int *n, const unsigned int len, const char *name ) {
	unsigned int i;
	double *y;
	FILE *f;

	// Allocate an output array:
	y = (double*) malloc( len * sizeof(double) );
	if ( y == NULL ) {
		printf( "Error allocating memory.\n" );
		exit( 1 );
	}

	// Generate fixture data:
	for ( i = 0; i < len; i++ ) {
		y[ i ] = boost::math::falling_factorial( x[ i ], n[ i ] );
	}
	// Open a new file:
	f = fopen( name, "w" );
	if ( f == NULL ) {
		printf( "Error opening file.\n" );
		exit( 1 );
	}
	// Write data as JSON:
	write_data_as_json( f, x, n, y, len );

	// Close the file:
	fclose( f );

	// Free allocated memory:
	free( y );
}

/**
* Main execution sequence.
*/
int main( void ) {
	unsigned int len;
	double *x;
	int *n;

	// Define the array length:
	len = 1000;

	// Allocate arrays:
	x = (double*) malloc( len * sizeof(double) );
	if ( x == NULL ) {
		printf( "Error allocating memory.\n" );
		exit( 1 );
	}
	n = (int*) malloc( len * sizeof(int) );
	if ( n == NULL ) {
		printf( "Error allocating memory.\n" );
		exit( 1 );
	}

	// Generate fixture data:
	rand_array_i32( n, len, 1, 100 );

	linspace_f64( x, len, 0.0, 0.5 );
	generate( x, n, len, "small.json" );

	linspace_f64( x, len, 1.0, 40.0 );
	generate( x, n, len, "large.json" );

	linspace_f64( x, len, -0.1, -40.0 );
	generate( x, n, len, "negative.json" );

	// Free allocated memory:
	free( x );
	free( n );

	return 0;
}

max_stars_repo_name: ghalimi/stdlib


text:
/*
  Copyright (c) 2014-2015 DataStax

  Licensed under the Apache License, Version 2.0 (the "License");
  you may not use this file except in compliance with the License.
  You may obtain a copy of the License at

  http://www.apache.org/licenses/LICENSE-2.0

  Unless required by applicable law or agreed to in writing, software
  distributed under the License is distributed on an "AS IS" BASIS,
  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  See the License for the specific language governing permissions and
  limitations under the License.
*/

#ifdef STAND_ALONE
#   define BOOST_TEST_MODULE cassandra
#endif

#include <algorithm>

#include <boost/test/unit_test.hpp>
#include <boost/test/debug.hpp>
#include <boost/lexical_cast.hpp>
#include <boost/cstdint.hpp>
#include <boost/format.hpp>
#include <boost/thread.hpp>
#include <boost/chrono.hpp>
#include <boost/atomic.hpp>
#include <boost/scoped_ptr.hpp>

#include "cassandra.h"
#include "test_utils.hpp"
#include "cql_ccm_bridge.hpp"

struct Version1DowngradeTests {
  Version1DowngradeTests() {}
};

BOOST_FIXTURE_TEST_SUITE(version1_downgrade, Version1DowngradeTests)

BOOST_AUTO_TEST_CASE(query_after_downgrade)
{
  test_utils::CassLog::reset("Error response: 'Invalid or unsupported protocol version: 2");

  size_t row_count;

  {
    test_utils::CassClusterPtr cluster(cass_cluster_new());

    const cql::cql_ccm_bridge_configuration_t& conf = cql::get_ccm_bridge_configuration("config_v1.txt");

    boost::shared_ptr<cql::cql_ccm_bridge_t> ccm = cql::cql_ccm_bridge_t::create_and_start(conf, "test", 3);

    test_utils::initialize_contact_points(cluster.get(), conf.ip_prefix(), 3, 0);

    cass_cluster_set_protocol_version(cluster.get(), 2);

    test_utils::CassSessionPtr session(test_utils::create_session(cluster.get()));

    test_utils::CassResultPtr result;
    test_utils::execute_query(session.get(), "SELECT * FROM system.schema_keyspaces", &result);

    row_count = cass_result_row_count(result.get());
  }

  BOOST_CHECK(row_count > 0);
  BOOST_CHECK(test_utils::CassLog::message_count() > 0);
}

BOOST_AUTO_TEST_SUITE_END()

max_stars_repo_name: stanislav-podlesny/cpp-driver

max_issues_repo_name: stanislav-podlesny/cpp-driver

max_forks_repo_name: stanislav-podlesny/cpp-driver


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 PositionsFile.cpp
 *
 * \brief Implementation of positions file.
 *
 * \author Marek Balint
 * \author Mauricio Varea
 * \author Aaron Day
 */

#include "data/PositionsFile.hh"

#include <vector>
#include <algorithm>
#include <numeric>
#include <utility>

#include <boost/format.hpp>
#include <boost/foreach.hpp>
#include <boost/filesystem.hpp>
#include <boost/filesystem/operations.hpp>
#include <boost/shared_ptr.hpp>
#include <boost/iostreams/stream.hpp>
#include <boost/iostreams/device/array.hpp>
#include <boost/algorithm/string.hpp>

namespace bcl2fastq {
namespace data {

PositionsFile::PositionsFile(const common::RawDataBuffer& data,
                             bool ignoreErrors /*= false*/
)
    : FileReaderBase(data,
                     ignoreErrors)
{
}

namespace detail {

///////////////////////////////////////////////////////////////////////////////
// positions file type: LOCS
///////////////////////////////////////////////////////////////////////////////

/// \brief LOCS file.
class LocsFile : public PositionsFile, public BinaryAllClustersFileReader
{

public:

    static const boost::filesystem::path EXTENSION;

    /// \brief Constructor.
    /// \param data Raw data from the LOCS file.
    /// \param ignoreErrors Suppress errors opening file and/or reading it.
    explicit LocsFile(
        const common::RawDataBuffer& data,
        bool ignoreErrors,
        bool skipHeader
    );

    /// \brief Read records from the file to buffer.
    /// \param targetBuffer Target buffer to read to.
    /// \param targetSize Maximum number of records to be read.
    /// \return Number of records read.
    virtual std::size_t read(
        std::vector<data::PositionsFile::Record>& targetBuffer,
        std::size_t targetSize
    );

private:

    /// \brief LOCS header type definition.
#pragma pack(push, 1)
    struct Header
    {
        /// \brief Positions file clusters count.
        typedef uint32_t PositionsFileClustersCount;

        /// \brief Value of the first constant field.
        static const uint32_t FIELD1;

        /// \brief Value of the second constant field.
        static const float FIELD2;

        /// \brief The first constant field.
        uint32_t field1_;

        /// \brief The second constant field.
        float field2_;

        /// \brief Clusters count.
        PositionsFileClustersCount clustersCount_;
    };

    /// \brief LOCS record type definition
    struct Record
    {
        /// \brief X-coordinate.
        float x_;

        /// \brief y-coordinate.
        float y_;

    };
#pragma pack(pop)

    /// \brief Get the number of bytes in a record
    /// \return Number of bytes in a record
    virtual std::size_t getRecordBytes() const { return sizeof(Record); }

    /// \brief Read the header
    virtual bool readHeader(bool ignoreErrors);

    /// \brief Validate the header. Throw an exception on failure.
    virtual bool validateHeader(Header& header, bool ignoreErrors);

    /// \brief Gets the file type string used for error messages.
    /// \return File type string
    virtual std::string getFileTypeStr() const { return "LOCS"; }

    /// \brief Transform the Record type into a PositionsFile::Record
    /// \param record Record to transform
    static PositionsFile::Record transformRecord(const Record &record);
};

const boost::filesystem::path LocsFile::EXTENSION(".locs");

LocsFile::LocsFile(
    const common::RawDataBuffer& data,
    bool ignoreErrors,
    bool skipHeader
)
: FileReaderBase(data,
                 ignoreErrors)
, PositionsFile(data,
                ignoreErrors)
, BinaryAllClustersFileReader(data,
                              ignoreErrors)
{
    if (!skipHeader)
    {
        readHeader(ignoreErrors);
    }
}

bool LocsFile::readHeader(bool ignoreErrors)
{
    Header header;
    if (!BinaryAllClustersFileReader::readHeader(header) ||
        !validateHeader(header, ignoreErrors))
    {
        return false;
    }

    clustersCount_ = static_cast<common::ClustersCount>(header.clustersCount_);

    return true;
}

bool LocsFile::validateHeader(LocsFile::Header& header, bool ignoreErrors)
{
    int errnum = errno;
    if (header.field1_ != Header::FIELD1)
    {
        const uint32_t field1 = header.field1_;
        const std::string message =
            (boost::format("Corrupted header of LOCS file '%s': header_field_1=%d")
                        % this->getPath().string() % field1).str();
        if (ignoreErrors) {
            BCL2FASTQ_LOG(common::LogLevel::WARNING) << message << std::endl;
        }
        else {
            BOOST_THROW_EXCEPTION(common::InputDataError(errnum, message));
        }

        clustersCount_ = 0;
        return false;
    }
    if (header.field2_ != Header::FIELD2)
    {
        const float field2 = header.field2_;
        const std::string message =
            (boost::format("Corrupted header of LOCS file '%s': header_field_2=%f")
                        % this->getPath().string() % field2).str();

        if (ignoreErrors) {
            BCL2FASTQ_LOG(common::LogLevel::WARNING) << message << std::endl;
        }
        else {
            BOOST_THROW_EXCEPTION(common::InputDataError(errnum, message));
        }

        clustersCount_ = 0;
        return false;
    }

    clustersCount_ = static_cast<common::ClustersCount>(header.clustersCount_);

    return true;
}

std::size_t LocsFile::read(
    std::vector<data::PositionsFile::Record>& targetBuffer,
    std::size_t targetSize
)
{
    targetBuffer.resize(targetSize);
    std::size_t ret = readBytes(reinterpret_cast<char*>(targetBuffer.data()), targetSize*sizeof(Record));

    const Record * const begin = reinterpret_cast<const Record * const>(targetBuffer.data());
    const Record * const end = begin + targetSize;

    std::transform(begin, end, targetBuffer.data(), &LocsFile::transformRecord);

    return ret;
}

const uint32_t LocsFile::Header::FIELD1 = 1;

const float LocsFile::Header::FIELD2 = 1.0;

PositionsFile::Record LocsFile::transformRecord(const Record &record)
{
    PositionsFile::Record ret;
    ret.x_ = static_cast<PositionsFile::Record::ClusterCoordinate>(record.x_ * 10.0 + 1000.5);
    ret.y_ = static_cast<PositionsFile::Record::ClusterCoordinate>(record.y_ * 10.0 + 1000.5);

    return ret;
}

///////////////////////////////////////////////////////////////////////////////
// positions file type: CLOCS
///////////////////////////////////////////////////////////////////////////////

/// \brief CLOCS file.
class ClocsFile : public PositionsFile, public BinaryAllClustersFileReader
{
public:

    static const boost::filesystem::path EXTENSION;

    /// \brief Constructor.
    /// \param data Raw data from the CLOCS file.
    /// \param ignoreErrors Suppress errors opening file and/or reading it.
    ClocsFile(
        const common::RawDataBuffer& data,
        bool ignoreErrors,
        bool skipHeader
    );

    /// \brief Read records from the file to buffer.
    /// \param targetBuffer Target buffer to read to.
    /// \param targetSize Maximum number of records to be read.
    /// \return Number of records read.
    virtual std::size_t read(
        std::vector<data::PositionsFile::Record>& targetBuffer,
        std::size_t targetSize
    );

private:

    /// \brief CLOCS header type definition.
#pragma pack(push, 1)
    struct Header
    {
    public:
        /// \brief supported version number
        static const uint8_t VERSION;

        /// \brief Version of the format.
        uint8_t version_;

        /// \brief Number of tiles
        uint32_t numTiles_;
    };

    struct Record
    {
    public:

        /// \brief x-offset
        uint8_t xOffset_;

        /// \brief y-offset
        uint8_t yOffset_;

    };
#pragma pack(pop)

private:

    /// \brief Get the number of bytes in a record
    /// \return Number of bytes in a record
    virtual std::size_t getRecordBytes() const { return sizeof(Record); }

    /// \brief Validate the header. Throw an exception on failure.
    virtual bool readHeader(bool ignoreErrors);

    /// \brief Validate the header. Throw an exception on failure.
    virtual bool validateHeader(Header& header, bool ignoreErrors);

//    /// \brief Read the number of clusters in the tile.
//    uint8_t readClusterCount();

    /// \brief Gets the file type string used for error messages.
    /// \return File type string
    virtual std::string getFileTypeStr() const { return "CLOCS"; }

    /// \brief Transform the Record type into a PositionsFile::Record
    /// \param record Record to transform
    PositionsFile::Record transformRecord(const Record &record);

private:

    /// \brief Number of tiles
    uint32_t numTiles_;

    uint32_t tileIndex_;
};

const boost::filesystem::path ClocsFile::EXTENSION(".clocs");

ClocsFile::ClocsFile(
    const common::RawDataBuffer& data,
    bool ignoreErrors,
    bool skipHeader
)
: FileReaderBase(data,
                 ignoreErrors)
, PositionsFile(data,
                ignoreErrors)
, BinaryAllClustersFileReader(data,
                              ignoreErrors)
, numTiles_(0)
, tileIndex_(0)
{
    if (!skipHeader)
    {
        readHeader(ignoreErrors);
    }
}

bool ClocsFile::readHeader(bool ignoreErrors)
{
    Header header;
    if (!BinaryAllClustersFileReader::readHeader(header) ||
        !validateHeader(header, ignoreErrors))
    {
        return false;
    }

    numTiles_ = header.numTiles_;

    return true;
}

bool ClocsFile::validateHeader(ClocsFile::Header& header,
                               bool ignoreErrors)
{
    int errnum = errno;
    if (header.version_ != Header::VERSION) {

        const std::string message =
            (boost::format("Corrupted header of CLOCS file '%s': header_version=%d")
                    % this->getPath().string() % header.version_).str();

        if (ignoreErrors) {
           BCL2FASTQ_LOG(common::LogLevel::WARNING) << message << std::endl;
        }
        else {
            BOOST_THROW_EXCEPTION(common::InputDataError(errnum, message));
        }

        return false;
    }

    return true;
}

//uint8_t ClocsFile::readClusterCount()
//{
//    uint8_t clustersCount = 0;
//
//    std::streamsize clustersCountLength =
//        BinaryFileReader::read(reinterpret_cast<char *>(&clustersCount),
//                               sizeof(clustersCount));
//
//    if (clustersCountLength != sizeof(clustersCount))
//    {
//        logError(clustersCountLength,
//                 sizeof(clustersCount));
//    }
//
//    return clustersCount;
//}

std::size_t ClocsFile::read(
    std::vector<data::PositionsFile::Record>& targetBuffer,
    std::size_t targetSize
)
{
    targetBuffer.clear();

    size_t numClustersRead = 0;
    if (!validateCondition(numTiles_ > 0, "No tiles read from clocs file.")) { return 0; }

    // Skip first 5 bytes, which are the version number and number of bins in the file
    size_t byteIndex = 5;
    for (uint32_t binNumber = 0; binNumber < numTiles_; ++binNumber)
    {
        if (!validateCondition(data_->size() > byteIndex, "Malformed clocs file.")) { return 0; }
        unsigned char numClusters = static_cast<unsigned char>(data_->at(byteIndex++));

        if (numClusters == 0)
        {
            continue;
        }

        targetBuffer.resize(targetBuffer.size()+numClusters);
        if (!validateCondition(data_->size() >= byteIndex + numClusters*2, "Malformed clocs file.")) { return 0; }
        for (uint32_t i = 0; i < numClusters; ++i)
        {
            auto& record = targetBuffer[numClustersRead + i];
            record.x_ = static_cast<PositionsFile::Record::ClusterCoordinate>( static_cast<uint8_t>((*data_)[byteIndex++]) + (binNumber % 82)*250) + 1000.5;
            record.y_ = static_cast<PositionsFile::Record::ClusterCoordinate>( static_cast<uint8_t>((*data_)[byteIndex++]) + (binNumber / 82)*250) + 1000.5;
        }

        numClustersRead += numClusters;
    }

    return numClustersRead;
}

const uint8_t ClocsFile::Header::VERSION = 0x01;

PositionsFile::Record ClocsFile::transformRecord(const Record &record)
{
    PositionsFile::Record ret;
    ret.x_ = static_cast<PositionsFile::Record::ClusterCoordinate>(record.xOffset_ + (tileIndex_ % 82)*250) + 1000.5;
    ret.y_ = static_cast<PositionsFile::Record::ClusterCoordinate>(record.yOffset_ + (tileIndex_ / 82)*250) + 1000.5;

    return ret;
}

class PosFile : public PositionsFile, FileReader
{
public:
    static const boost::filesystem::path EXTENSION;

    PosFile(const common::RawDataBuffer& data,
            bool                         ignoreErrors);

    /// \brief Read records from file to buffer.
    /// \param targetBuffer Target buffer to read to.
    /// \param targetSize Maximum number of records to be read.
    /// \return Number of records read.
    virtual std::size_t read(
        std::vector<data::PositionsFile::Record>& targetBuffer,
        std::size_t targetSize
    );

private:
    /// \brief Gets the file type string used for error messages.
    /// \return File type string
    virtual std::string getFileTypeStr() const { return "POS"; }
};

const boost::filesystem::path PosFile::EXTENSION(".txt");

PosFile::PosFile(const common::RawDataBuffer& data,
                 bool                         ignoreErrors)
: FileReaderBase(data,
                 ignoreErrors)
, PositionsFile(data,
                ignoreErrors)
, FileReader(data,
             ignoreErrors)
{
}

std::size_t PosFile::read(
    std::vector<data::PositionsFile::Record>& targetBuffer,
    std::size_t targetSize
)
{
    targetBuffer.clear();
    targetBuffer.reserve(targetSize);

    std::size_t recordsRead = 0;
    static const std::string errMsgBegin = "Ignoring read failure for POS file '";
    try
    {
        std::string line;
        while (std::getline(istr_, line) && recordsRead < targetSize)
        {
            std::vector<std::string> strs;
            boost::split(strs, line, boost::is_any_of("\t "));

            if (strs.size() != 2)
            {
                std::string errMsg = "Formatting error for POS file: " + this->getPath().string();
                if (ignoreErrors_)
                {
                    BCL2FASTQ_LOG(common::LogLevel::WARNING) << errMsg << std::endl;
                    return recordsRead;
                }
                else
                {
                    int errnum = errno;
                    BOOST_THROW_EXCEPTION(common::InputDataError(errnum, errMsg));
                }
            }

            PositionsFile::Record record;
            record.x_ = boost::lexical_cast<PositionsFile::Record::ClusterCoordinate>(strs[0]);
            record.y_ = boost::lexical_cast<PositionsFile::Record::ClusterCoordinate>(strs[1]);
            targetBuffer.push_back(record);

            ++recordsRead;
        }
    }
    CATCH_AND_IGNORE(common::IoError, errMsgBegin)
    CATCH_AND_IGNORE(std::ios_base::failure, errMsgBegin)
    CATCH_AND_IGNORE_ALL(errMsgBegin)

    return recordsRead;
}

} // namespace detail

bool PositionsFileFactory::doesFileExist(
    const boost::filesystem::path& intensitiesDir,
    bool                           aggregateTilesFlag,
    bool                           isPatternedFlowcell,
    common::LaneNumber             laneNumber,
    common::TileNumber             tileNumber,
    size_t&                        headerSize,
    boost::filesystem::path&       positionsFilePath
)
{
    if (binaryFilePathExists(intensitiesDir,
                             aggregateTilesFlag,
                             isPatternedFlowcell,
                             laneNumber,
                             tileNumber,
                             detail::LocsFile::EXTENSION,
                             positionsFilePath))
    {
        headerSize = 12;
        return true;
    }

    if (binaryFilePathExists(intensitiesDir,
                             aggregateTilesFlag,
                             isPatternedFlowcell,
                             laneNumber,
                             tileNumber,
                             detail::ClocsFile::EXTENSION,
                             positionsFilePath))
    {
        headerSize = 5;
        return true;
    }

    if (posFilePathExists(intensitiesDir,
                          laneNumber,
                          tileNumber,
                          positionsFilePath))
    {
        headerSize = 0;
        return true;
    }

    return false;
}

bool PositionsFileFactory::binaryFilePathExists(
    const boost::filesystem::path& intensitiesDir,
    bool                           aggregateTilesFlag,
    bool                           isPatternedFlowcell,
    common::LaneNumber             laneNumber,
    common::TileNumber             tileNumber,
    const boost::filesystem::path& fileExtension,
    boost::filesystem::path&       positionsFilePath
)
{
    if( isPatternedFlowcell )
    {
        positionsFilePath = boost::filesystem::path(
            intensitiesDir
            /
            boost::filesystem::path("s" + fileExtension.string())
        );
    } else if( aggregateTilesFlag ) {
        positionsFilePath = boost::filesystem::path(
            intensitiesDir
            /
            boost::filesystem::path((boost::format("L%03d") % laneNumber).str())
            /
            boost::filesystem::path((boost::format("s_%d" + fileExtension.string()) % laneNumber).str())
        );
    } else {
        positionsFilePath = boost::filesystem::path(
            intensitiesDir
            /
            boost::filesystem::path((boost::format("L%03d") % laneNumber).str())
            /
            boost::filesystem::path((boost::format("s_%d_%d" + fileExtension.string()) % laneNumber % tileNumber).str())
        );
    }

    return boost::filesystem::exists(positionsFilePath);
}

bool PositionsFileFactory::posFilePathExists(
    const boost::filesystem::path& intensitiesDir,
    common::LaneNumber             laneNumber,
    common::TileNumber             tileNumber,
    boost::filesystem::path&       positionsFilePath
)
{
    positionsFilePath = boost::filesystem::path(
        intensitiesDir
        /
        boost::filesystem::path((boost::format("s_%d_%03d_pos.txt") % laneNumber % tileNumber).str())
    );

    return boost::filesystem::exists(positionsFilePath);
}

std::streamsize PositionsFileFactory::read(const common::RawDataBuffer& inputData,
                                           bool ignoreErrors,
                                           bool skipHeader,
                                           std::vector<data::PositionsFile::Record>& targetBuffer,
                                           std::size_t targetSize)
{
    auto extension = inputData.path_.extension();
    if (extension == detail::LocsFile::EXTENSION)
    {
        detail::LocsFile posFile(inputData,
                                 ignoreErrors,
                                 skipHeader);

        return posFile.read(targetBuffer, targetSize > 0 ? targetSize : posFile.getClustersCount());

    }
    else if (extension == detail::ClocsFile::EXTENSION)
    {
        detail::ClocsFile posFile(inputData,
                                  ignoreErrors,
                                  skipHeader);

        return posFile.read(targetBuffer, targetSize);
    }
    else if (extension == detail::PosFile::EXTENSION)
    {
        detail::PosFile posFile(inputData,
                                ignoreErrors);

        posFile.read(targetBuffer, targetSize);
    }

    return 0;
}

} // namespace data
} // namespace bcl2fastq


text:
// ----------------------------------------------------------------------------
// -                        Open3D: www.open3d.org                            -
// ----------------------------------------------------------------------------
// The MIT License (MIT)
//
// Copyright (c) 2018 www.open3d.org
//
// Permission is hereby granted, free of charge, to any person obtaining a copy
// of this software and associated documentation files (the "Software"), to deal
// in the Software without restriction, including without limitation the rights
// to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the Software is
// furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included in
// all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
// IN THE SOFTWARE.
// ----------------------------------------------------------------------------

#include "open3d/t/geometry/TriangleMesh.h"

#include <Eigen/Core>
#include <string>
#include <unordered_map>

#include "open3d/core/EigenConverter.h"
#include "open3d/core/ShapeUtil.h"
#include "open3d/core/Tensor.h"
#include "open3d/core/TensorList.h"

namespace open3d {
namespace t {
namespace geometry {

TriangleMesh::TriangleMesh(core::Dtype vertex_dtype,
                           core::Dtype triangle_dtype,
                           const core::Device &device)
    : Geometry(Geometry::GeometryType::TriangleMesh, 3),
      device_(device),
      vertex_attr_(TensorListMap("vertices")),
      triangle_attr_(TensorListMap("triangles")) {
    SetVertices(core::TensorList({3}, vertex_dtype, device_));
    SetTriangles(core::TensorList({3}, triangle_dtype, device_));
}

TriangleMesh::TriangleMesh(const core::TensorList &vertices,
                           const core::TensorList &triangles)
    : TriangleMesh(vertices.GetDtype(), triangles.GetDtype(), [&]() {
          if (vertices.GetDevice() != triangles.GetDevice()) {
              utility::LogError(
                      "vertices' device {} does not match triangles' device "
                      "{}.",
                      vertices.GetDevice().ToString(),
                      triangles.GetDevice().ToString());
          }
          return vertices.GetDevice();
      }()) {
    SetVertices(vertices);
    SetTriangles(triangles);
}

}  // namespace geometry
}  // namespace t
}  // namespace open3d