Gradient.h

Go to the documentation of this file.

#ifndef GRADIENT_H_
#define GRADIENT_H_

#include <string>
#include "Vector.h"
#include "Data.h"

class Gradient {
public:
        Gradient(Data* input, Vector<> boundaryValues = Vector<>(0,0,0));
        ~Gradient();

        // Sizes in x/y/z dimension
        unsigned short xdim() const { return mX; }
        unsigned short ydim() const { return mY; }
        unsigned short zdim() const { return mZ; }

        Vector<> getData(unsigned int aX, unsigned int aY, unsigned int aZ) const{
                assert(aX < mX);
                assert(aY < mY);
                assert(aZ < mZ);

                const unsigned int Ystride = mX;
                const unsigned int Zstride = Ystride * mY;

                return mData[aZ * Zstride + aY * Ystride + aX];
        }
        
        Vector<> getData_RangeCheck(int aX, int aY, int aZ){
                
        //      printf("%ix%ix%i\n",aX,aY,aZ);
                if ((aX<0) || (aX>=mX)) return Vector<>(0,0,0);
                if ((aY<0) || (aY>=mY)) return Vector<>(0,0,0);
                if ((aZ<0) || (aZ>=mZ)) return Vector<>(0,0,0);

//              return getData((unsigned int) aX, (unsigned int) aY, (unsigned int) aZ);
                const unsigned int Ystride = (unsigned int) mX;
                const unsigned int Zstride = (unsigned int) Ystride * mY;

                return mData[aZ * Zstride + aY * Ystride + aX];         
        }
        
        Vector<> getDataTrilinearInterpolated(float aX, float aY, float aZ){
/*
Trilinear interpolation. cube-corners:
  E-----F
 /|    /| 
A-----B |                               Y  Z
| H---|-G                               | /
|/        |/                            |/
D-----C                                 ---X

1.Step: 4 * interpolate Z-Pos.
2.Step: 2 * interpolate Y-Pos.
3.Step: 1 *     interpolate X-Pos.

*/
                        
        //      printf("Input: %f %f %f\n",aX,aY,aZ);
                int lowerX = (int) floor(aX);
                int lowerY = (int) floor(aY);
                int lowerZ = (int) floor(aZ);
                int upperX = lowerX+1;
                int upperY = lowerY+1;
                int upperZ = lowerZ+1;
        //      printf("Output: %i %i %i\n", lowerX,lowerY,lowerZ);
        //define cube corners:
                Vector<> A = getData_RangeCheck(lowerX, upperY, lowerZ);
                Vector<> B = getData_RangeCheck(upperX, upperY, lowerZ);                
                Vector<> C = getData_RangeCheck(upperX, lowerY, lowerZ);
                Vector<> D = getData_RangeCheck(lowerX, lowerY, lowerZ);

                Vector<> E = getData_RangeCheck(lowerX, upperY, upperZ);
                Vector<> F = getData_RangeCheck(upperX, upperY, upperZ);
                Vector<> G = getData_RangeCheck(upperX, lowerY, upperZ);
                Vector<> H = getData_RangeCheck(lowerX, lowerY, upperZ);
        
                float floatX = aX - lowerX;
                float floatY = aY - lowerY;
                float floatZ = aZ - lowerZ;
        //      printf("Floats: %f %f %f\n\n",floatX,floatY,floatZ);            

//      1.Step:
        //      Vector<> AE = Vector<> (E.x() * floatZ + A.x() * (1-floatZ) , E.y() * floatZ + A.y() * (1-floatZ) , E.z() * floatZ + A.z() * (1-floatZ) );
        //      Vector<> BF = Vector<> (F.x() * floatZ + B.x() * (1-floatZ) , F.y() * floatZ + B.y() * (1-floatZ) , F.z() * floatZ + B.z() * (1-floatZ) );
        //      Vector<> CG = Vector<> (G.x() * floatZ + C.x() * (1-floatZ) , G.y() * floatZ + C.y() * (1-floatZ) , G.z() * floatZ + C.z() * (1-floatZ) );

                Vector<> AE = E * floatZ + A * (1-floatZ);
                Vector<> BF = F * floatZ + B * (1-floatZ);
                Vector<> CG = G * floatZ + C * (1-floatZ);

                Vector<> DH = H * floatZ + D * (1-floatZ);
//      2.Step:
                Vector<> AEDH = AE * floatY + DH * (1-floatY);
                Vector<> BFCG = BF * floatY + CG * (1-floatY);
//      3.Step:
                return (BFCG * floatX + AEDH * (1-floatX)).normalized();
        }

        
        

private:
        Vector<>* mData;
        
        Vector<> mBoundaryValues;

        unsigned short mX, mY, mZ;
        
};

#endif

---