ThreeDGrid.h

Go to the documentation of this file.

/*
 *  Copyright (c) 2012, Chris Rogers
 *  All rights reserved.
 *  Redistribution and use in source and binary forms, with or without
 *  modification, are permitted provided that the following conditions are met:
 *  1. Redistributions of source code must retain the above copyright notice,
 *     this list of conditions and the following disclaimer.
 *  2. 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.
 *  3. Neither the name of STFC 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 OWNER 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.
 */
 
#ifndef _OPALX_FIELDS_THREEDGRID_HH_
#define _OPALX_FIELDS_THREEDGRID_HH_
 
#include <algorithm>
#include <cmath>
#include <vector>
 
#include "Fields/Interpolation/Mesh.h"
 
namespace interpolation {
 
    class ThreeDGrid : public Mesh {
    public:
        class Iterator;
 
        ThreeDGrid* clone() { return new ThreeDGrid(*this); }
 
        Mesh* dual() const;
 
        ThreeDGrid();
 
        ThreeDGrid(
                double dX, double dY, double dZ, double minX, double minY, double minZ,
                int numberOfXCoords, int numberOfYCoords, int numberOfZCoords);
 
        ThreeDGrid(
                int xSize, const double* x, int ySize, const double* y, int zSize, const double* z);
 
        ThreeDGrid(std::vector<double> x, std::vector<double> y, std::vector<double> z);
 
        ~ThreeDGrid();
 
        inline double& x(const int& i) { return x_m[i - 1]; }
 
        inline double& y(const int& j) { return y_m[j - 1]; }
 
        inline double& z(const int& k) { return z_m[k - 1]; }
 
        inline const double& x(const int& i) const { return x_m[i - 1]; }
 
        inline const double& y(const int& j) const { return y_m[j - 1]; }
 
        inline const double& z(const int& j) const { return z_m[j - 1]; }
 
        inline int xSize() const { return static_cast<int>(x_m.size()); }
 
        inline int ySize() const { return static_cast<int>(y_m.size()); }
 
        inline int zSize() const { return static_cast<int>(z_m.size()); }
 
        std::vector<double> xVector() { return std::vector<double>(x_m); }
 
        std::vector<double> yVector() { return std::vector<double>(y_m); }
 
        std::vector<double> zVector() { return std::vector<double>(z_m); }
 
        inline double* newXArray();
 
        inline double* newYArray();
 
        inline double* newZArray();
 
        inline void xLowerBound(const double& x, int& xIndex) const;
 
        inline void yLowerBound(const double& y, int& yIndex) const;
 
        inline void zLowerBound(const double& z, int& zIndex) const;
 
        inline void lowerBound(
                const double& x, int& xIndex, const double& y, int& yIndex, const double& z,
                int& zIndex) const;
 
        inline void lowerBound(
                const double& x, const double& y, const double& z, Mesh::Iterator& it) const;
 
        inline double minX() const;
 
        inline double maxX() const;
 
        inline double minY() const;
 
        inline double maxY() const;
 
        inline double minZ() const;
 
        inline double maxZ() const;
 
        inline void setX(int nXCoords, double* x);
 
        inline void setY(int nYCoords, double* y);
 
        inline void setZ(int nZCoords, double* z);
 
        void add(VectorMap* map);
 
        void remove(VectorMap* map);
 
        Mesh::Iterator begin() const;
 
        Mesh::Iterator end() const;
 
        virtual void getPosition(const Mesh::Iterator& it, double* position) const;
 
        inline int getPositionDimension() const;
 
        inline int toInteger(const Mesh::Iterator& lhs) const;
 
        inline void setConstantSpacing(bool spacing);
 
        void setConstantSpacing();
 
        inline bool getConstantSpacing() const;
 
        Mesh::Iterator getNearest(const double* position) const;
 
        static void vectorLowerBound(std::vector<double> vec, double x, int& index);
 
    protected:
        // Change position
        virtual Mesh::Iterator& addEquals(Mesh::Iterator& lhs, int difference) const;
        virtual Mesh::Iterator& subEquals(Mesh::Iterator& lhs, int difference) const;
        virtual Mesh::Iterator& addEquals(Mesh::Iterator& lhs, const Mesh::Iterator& rhs) const;
        virtual Mesh::Iterator& subEquals(Mesh::Iterator& lhs, const Mesh::Iterator& rhs) const;
        virtual Mesh::Iterator& addOne(Mesh::Iterator& lhs) const;
        virtual Mesh::Iterator& subOne(Mesh::Iterator& lhs) const;
        // Check relative position
        virtual bool isGreater(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs) const;
 
    private:
        std::vector<double> x_m;
        std::vector<double> y_m;
        std::vector<double> z_m;
        int xSize_m;
        int ySize_m;
        int zSize_m;
        std::vector<VectorMap*> maps_m;
        bool constantSpacing_m;
 
        friend Mesh::Iterator operator++(Mesh::Iterator& lhs, int);
        friend Mesh::Iterator operator--(Mesh::Iterator& lhs, int);
        friend Mesh::Iterator& operator++(Mesh::Iterator& lhs);
        friend Mesh::Iterator& operator--(Mesh::Iterator& lhs);
        friend Mesh::Iterator operator-(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
        friend Mesh::Iterator operator+(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
        friend Mesh::Iterator& operator-=(Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
        friend Mesh::Iterator& operator+=(Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
 
        friend bool operator==(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
        friend bool operator!=(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
        friend bool operator>=(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
        friend bool operator<=(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
        friend bool operator<(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
        friend bool operator>(const Mesh::Iterator& lhs, const Mesh::Iterator& rhs);
    };
 
    double* ThreeDGrid::newXArray() {
        double* x = new double[x_m.size()];
        for (unsigned int i = 0; i < x_m.size(); i++)
            x[i] = x_m[i];
        return x;
    }
 
    double* ThreeDGrid::newYArray() {
        double* y = new double[y_m.size()];
        for (unsigned int i = 0; i < y_m.size(); i++)
            y[i] = y_m[i];
        return y;
    }
 
    double* ThreeDGrid::newZArray() {
        double* z = new double[z_m.size()];
        for (unsigned int i = 0; i < z_m.size(); i++)
            z[i] = z_m[i];
        return z;
    }
 
    void ThreeDGrid::lowerBound(
            const double& x, int& xIndex, const double& y, int& yIndex, const double& z,
            int& zIndex) const {
        xLowerBound(x, xIndex);
        yLowerBound(y, yIndex);
        zLowerBound(z, zIndex);
    }
 
    void ThreeDGrid::lowerBound(
            const double& x, const double& y, const double& z, Mesh::Iterator& it) const {
        xLowerBound(x, it[0]);
        yLowerBound(y, it[1]);
        zLowerBound(z, it[2]);
        it[0]++;
        it[1]++;
        it[2]++;
    }
 
    void ThreeDGrid::xLowerBound(const double& x, int& xIndex) const {
        if (constantSpacing_m)
            xIndex = static_cast<int>(std::floor((x - x_m[0]) / (x_m[1] - x_m[0])));
        else {
            vectorLowerBound(x_m, x, xIndex);
        }
    }
 
    void ThreeDGrid::yLowerBound(const double& y, int& yIndex) const {
        if (constantSpacing_m)
            yIndex = static_cast<int>(std::floor((y - y_m[0]) / (y_m[1] - y_m[0])));
        else
            vectorLowerBound(y_m, y, yIndex);
    }
 
    void ThreeDGrid::zLowerBound(const double& z, int& zIndex) const {
        if (constantSpacing_m)
            zIndex = static_cast<int>(std::floor((z - z_m[0]) / (z_m[1] - z_m[0])));
        else
            vectorLowerBound(z_m, z, zIndex);
    }
 
    double ThreeDGrid::minX() const { return x_m[0]; }
 
    double ThreeDGrid::maxX() const { return x_m[xSize_m - 1]; }
 
    double ThreeDGrid::minY() const { return y_m[0]; }
 
    double ThreeDGrid::maxY() const { return y_m[ySize_m - 1]; }
 
    double ThreeDGrid::minZ() const { return z_m[0]; }
 
    double ThreeDGrid::maxZ() const { return z_m[zSize_m - 1]; }
 
    void ThreeDGrid::setX(int nXCoords, double* x) { x_m = std::vector<double>(x, x + nXCoords); }
 
    void ThreeDGrid::setY(int nYCoords, double* y) { y_m = std::vector<double>(y, y + nYCoords); }
 
    void ThreeDGrid::setZ(int nZCoords, double* z) { z_m = std::vector<double>(z, z + nZCoords); }
 
    int ThreeDGrid::getPositionDimension() const { return 3; }
 
    int ThreeDGrid::toInteger(const Mesh::Iterator& lhs) const {
        return (lhs.getState()[0] - 1) * zSize_m * ySize_m + (lhs.getState()[1] - 1) * zSize_m
               + (lhs.getState()[2] - 1);
    }
 
    void ThreeDGrid::setConstantSpacing(bool spacing) { constantSpacing_m = spacing; }
 
    bool ThreeDGrid::getConstantSpacing() const { return constantSpacing_m; }
}  // namespace interpolation
#endif