Commit 7f7a1237 authored by gsell's avatar gsell

bugfix in reading H5Block formatted field-maps and major code review

most code is now in a base class, all three classes are now using
the same interpolation code.
parent 3b011614
......@@ -11,6 +11,7 @@ set (_SRCS
StaticElectricField.cpp
StaticMagneticField.cpp
Fieldmap.cpp
FM3DH5BlockBase.cpp
FM3DH5Block.cpp
FM3DH5Block_nonscale.cpp
FM3DMagnetoStaticH5Block.cpp
......@@ -81,6 +82,7 @@ set (HDRS
FM2DElectroStatic.h
FM2DMagnetoStatic.h
FM3DDynamic.h
FM3DH5BlockBase.h
FM3DH5Block.h
FM3DH5Block_nonscale.h
FM3DMagnetoStaticH5Block.h
......
This diff is collapsed.
#ifndef CLASSIC_FIELDMAP3DH5BLOCK_HH
#define CLASSIC_FIELDMAP3DH5BLOCK_HH
//
// Class for dynamic 3D field-maps stored in H5hut files.
//
// Copyright (c) 2020, Achim Gsell, Paul Scherrer Institut, Villigen PSI, Switzerland
// All rights reserved.
//
// This file is part of OPAL.
//
// OPAL is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// You should have received a copy of the GNU General Public License
// along with OPAL. If not, see <https://www.gnu.org/licenses/>.
//
#ifndef CLASSIC_FIELDMAP3DH5BLOCK_H
#define CLASSIC_FIELDMAP3DH5BLOCK_H
#include "Fields/FM3DH5BlockBase.h"
#include "Fields/Fieldmap.h"
#include "hdf5.h"
#include "H5hut.h"
#include <vector>
class FM3DH5Block: public Fieldmap {
class FM3DH5Block: public FM3DH5BlockBase {
public:
virtual bool getFieldstrength(const Vector_t &R, Vector_t &E, Vector_t &B) const;
virtual void getFieldDimensions(double &zBegin, double &zEnd, double &rBegin, double &rEnd) const;
virtual void getFieldDimensions(double &xIni, double &xFinal, double &yIni, double &yFinal, double &zIni, double &zFinal) const;
virtual bool getFieldDerivative(const Vector_t &R, Vector_t &E, Vector_t &B, const DiffDirection &dir) const;
virtual void swap();
virtual void getInfo(Inform *msg);
virtual double getFrequency() const;
virtual void setFrequency(double freq);
virtual void getOnaxisEz(std::vector<std::pair<double, double> > & F);
virtual bool isInside(const Vector_t &r) const;
virtual bool getFieldstrength (
const Vector_t &R, Vector_t &E, Vector_t &B) const;
private:
FM3DH5Block(std::string aFilename);
~FM3DH5Block();
virtual void readMap();
virtual void freeMap();
std::vector<h5_float64_t> FieldstrengthEz_m; /**< 3D array with Ez */
std::vector<h5_float64_t> FieldstrengthEx_m; /**< 3D array with Ex */
std::vector<h5_float64_t> FieldstrengthEy_m; /**< 3D array with Ey */
std::vector<h5_float64_t> FieldstrengthHz_m; /**< 3D array with Hz */
std::vector<h5_float64_t> FieldstrengthHx_m; /**< 3D array with Hx */
std::vector<h5_float64_t> FieldstrengthHy_m; /**< 3D array with Hy */
h5_float64_t frequency_m;
FM3DH5Block (
std::string aFilename);
h5_float64_t xbegin_m;
h5_float64_t xend_m;
// int xcentral_idx_m;
virtual ~FM3DH5Block (
);
h5_float64_t ybegin_m;
h5_float64_t yend_m;
// int ycentral_idx_m;
virtual void readMap (
);
h5_float64_t zbegin_m;
h5_float64_t zend_m;
virtual void freeMap (
);
h5_float64_t hx_m; /**< length between points in grid, x-direction */
h5_float64_t hy_m; /**< length between points in grid, y-direction */
h5_float64_t hz_m; /**< length between points in grid, z-direction */
h5_int64_t num_gridpx_m; /**< Read in number of points after 0(not counted here) in grid, r-direction*/
h5_int64_t num_gridpy_m; /**< Read in number of points after 0(not counted here) in grid, r-direction*/
h5_int64_t num_gridpz_m; /**< Read in number of points after 0(not counted here) in grid, z-direction*/
std::vector<double> FieldstrengthHz_m; /**< 3D array with Hz */
std::vector<double> FieldstrengthHx_m; /**< 3D array with Hx */
std::vector<double> FieldstrengthHy_m; /**< 3D array with Hy */
bool swap_m;
friend class Fieldmap;
friend class FM3DH5BlockBase;
};
inline bool FM3DH5Block::isInside(const Vector_t &r) const
{
return ((r(0) >= xbegin_m && r(0) < xend_m) &&
(r(1) >= ybegin_m && r(1) < yend_m) &&
(r(2) >= zbegin_m && r(2) < zend_m));
}
#endif
#endif
\ No newline at end of file
// vi: set et ts=4 sw=4 sts=4:
// Local Variables:
// mode:c++
// c-basic-offset: 4
// indent-tabs-mode: nil
// require-final-newline: nil
// End:
This diff is collapsed.
//
// Base class for 3D field-maps in stored in H5hut files.
//
// Copyright (c) 2020, Achim Gsell, Paul Scherrer Institut, Villigen PSI, Switzerland
// All rights reserved.
//
// This file is part of OPAL.
//
// OPAL is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// You should have received a copy of the GNU General Public License
// along with OPAL. If not, see <https://www.gnu.org/licenses/>.
//
#ifndef CLASSIC_FIELDMAP3DH5BLOCKBASE_H
#define CLASSIC_FIELDMAP3DH5BLOCKBASE_H
#include "Fields/Fieldmap.h"
#include <vector>
#include "H5hut.h"
static_assert (sizeof(double) == sizeof (h5_float64_t),
"double and h5_float64_t are not the same type" );
static_assert (sizeof(long long) == sizeof (h5_int64_t),
"long long and h5_intt64_t are not the same type" );
class FM3DH5BlockBase: virtual public Fieldmap {
public:
virtual bool getFieldstrength (
const Vector_t& /*R*/, Vector_t& /*E*/, Vector_t& /*B*/) const {return false;};
void getFieldDimensions (
double &zBegin, double &zEnd,
double &rBegin, double &rEnd
) const {
zBegin = zbegin_m;
zEnd = zend_m;
rBegin = xbegin_m;
rEnd = xend_m;
}
void getFieldDimensions (
double &xIni, double &xFinal,
double &yIni, double &yFinal,
double &zIni, double &zFinal
) const {
xIni = xbegin_m;
xFinal = xend_m;
yIni = ybegin_m;
yFinal = yend_m;
zIni = zbegin_m;
zFinal = zend_m;
}
virtual bool getFieldDerivative (
const Vector_t &/*R*/,
Vector_t &/*E*/,
Vector_t &/*B*/,
const DiffDirection &/*dir*/
) const {
return false;
}
virtual void swap(
) {};
virtual void getInfo (
Inform *msg);
virtual double getFrequency (
) const;
virtual void setFrequency (
double freq);
void getOnaxisEz (
std::vector<std::pair<double, double> >& F);
protected:
FM3DH5BlockBase (
std::string aFilename) : Fieldmap (aFilename) {};
virtual ~FM3DH5BlockBase (
) {};
void OpenFileMPIOCollective (
const std::string aFilename);
long long GetNumSteps (
void);
void SetStep (
const long long);
void GetFieldInfo (
const char*);
void GetResonanceFrequency (
void);
void ReadField (
const char* name,
double* x,
double* y,
double* z
);
void CloseFile (
void);
virtual void readMap (
) {};
virtual void freeMap (
) {};
bool isInside (
const Vector_t &r
) const {
return ((r(0) >= xbegin_m && r(0) < xend_m) &&
(r(1) >= ybegin_m && r(1) < yend_m) &&
(r(2) >= zbegin_m && r(2) < zend_m));
}
struct IndexTriplet {
unsigned int i;
unsigned int j;
unsigned int k;
Vector_t weight;
IndexTriplet():
i(0),
j(0),
k(0),
weight(0.0)
{}
};
unsigned long getIndex (
unsigned int i,
unsigned int j,
unsigned int k
) const {
unsigned long result = i + j * num_gridpx_m;
result = k + result * num_gridpz_m;
return result;
}
IndexTriplet getIndex(const Vector_t &X) const {
IndexTriplet idx;
idx.i = std::floor((X(0) - xbegin_m) / hx_m);
idx.j = std::floor((X(1) - ybegin_m) / hy_m);
idx.k = std::floor((X(2) - zbegin_m) / hz_m);
PAssert_LT(idx.i, num_gridpx_m - 1);
PAssert_LT(idx.j, num_gridpy_m - 1);
PAssert_LT(idx.k, num_gridpz_m - 1);
idx.weight(0) = (X(0) - xbegin_m) / hx_m - idx.i;
idx.weight(1) = (X(1) - ybegin_m) / hy_m - idx.j;
idx.weight(2) = (X(2) - zbegin_m) / hz_m - idx.k;
return idx;
}
double getWeightedData (
const std::vector<double>& data,
const IndexTriplet& idx,
unsigned short corner) const;
Vector_t interpolateTrilinearly (
const std::vector<double>&,
const std::vector<double>&,
const std::vector<double>&,
const Vector_t& X) const;
enum { LX = 0, // low X
LY = 0, // low Y
LZ = 0, // low Z
HX = 4, // high X
HY = 2, // high Y
HZ = 1}; // high Z
h5_file_t file_m;
std::vector<double> FieldstrengthEz_m; /**< 3D array with Ez */
std::vector<double> FieldstrengthEx_m; /**< 3D array with Ex */
std::vector<double> FieldstrengthEy_m; /**< 3D array with Ey */
double xbegin_m;
double xend_m;
double ybegin_m;
double yend_m;
double zbegin_m;
double zend_m;
double hx_m; /**< length between points in grid, x-direction */
double hy_m; /**< length between points in grid, y-direction */
double hz_m; /**< length between points in grid, z-direction */
double num_gridpx_m; /**< number of points after 0(not counted here) in grid, r-direction*/
double num_gridpy_m; /**< number of points after 0(not counted here) in grid, r-direction*/
double num_gridpz_m; /**< number of points after 0(not counted here) in grid, z-direction*/
double frequency_m;
bool swap_m;
friend class Fieldmap;
};
#endif
// vi: set et ts=4 sw=4 sts=4:
// Local Variables:
// mode:c
// c-basic-offset: 4
// indent-tabs-mode: nil
// require-final-newline: nil
// End:
This diff is collapsed.
#ifndef CLASSIC_FIELDMAP3DH5BLOCK_nonscale_HH
#define CLASSIC_FIELDMAP3DH5BLOCK_nonscale_HH
#include "Fields/Fieldmap.h"
#include "hdf5.h"
#include "H5hut.h"
class FM3DH5Block_nonscale: public Fieldmap {
//
// Class for dynamic non-scaled 3D field-maps stored in H5hut files.
//
// Copyright (c) 2020, Achim Gsell, Paul Scherrer Institut, Villigen PSI, Switzerland
// All rights reserved.
//
// This file is part of OPAL.
//
// OPAL is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// You should have received a copy of the GNU General Public License
// along with OPAL. If not, see <https://www.gnu.org/licenses/>.
//
#ifndef CLASSIC_FIELDMAP3DH5BLOCK_NONSCALE_H
#define CLASSIC_FIELDMAP3DH5BLOCK_NONSCALE_H
#include "Fields/FM3DH5BlockBase.h"
#include <vector>
class FM3DH5Block_nonscale: public FM3DH5BlockBase {
public:
virtual bool getFieldstrength(const Vector_t &R, Vector_t &E, Vector_t &B) const;
virtual void getFieldDimensions(double &zBegin, double &zEnd, double &rBegin, double &rEnd) const;
virtual void getFieldDimensions(double &xIni, double &xFinal, double &yIni, double &yFinal, double &zIni, double &zFinal) const;
virtual bool getFieldDerivative(const Vector_t &R, Vector_t &E, Vector_t &B, const DiffDirection &dir) const;
virtual void swap();
virtual void getInfo(Inform *msg);
virtual double getFrequency() const;
virtual void setFrequency(double freq);
virtual void getOnaxisEz(std::vector<std::pair<double, double> > & F);
virtual bool isInside(const Vector_t &r) const;
virtual bool getFieldstrength (
const Vector_t &R, Vector_t &E, Vector_t &B) const;
private:
FM3DH5Block_nonscale(std::string aFilename);
~FM3DH5Block_nonscale();
virtual void readMap();
virtual void freeMap();
std::vector<h5_float64_t> FieldstrengthEz_m; /**< 3D array with Ez */
std::vector<h5_float64_t> FieldstrengthEx_m; /**< 3D array with Ex */
std::vector<h5_float64_t> FieldstrengthEy_m; /**< 3D array with Ey */
std::vector<h5_float64_t> FieldstrengthHz_m; /**< 3D array with Hz */
std::vector<h5_float64_t> FieldstrengthHx_m; /**< 3D array with Hx */
std::vector<h5_float64_t> FieldstrengthHy_m; /**< 3D array with Hy */
h5_float64_t frequency_m;
FM3DH5Block_nonscale (
std::string aFilename);
h5_float64_t xbegin_m;
h5_float64_t xend_m;
// int xcentral_idx_m;
virtual ~FM3DH5Block_nonscale (
);
h5_float64_t ybegin_m;
h5_float64_t yend_m;
// int ycentral_idx_m;
virtual void readMap (
);
h5_float64_t zbegin_m;
h5_float64_t zend_m;
virtual void freeMap (
);
h5_float64_t hx_m; /**< length between points in grid, x-direction */
h5_float64_t hy_m; /**< length between points in grid, y-direction */
h5_float64_t hz_m; /**< length between points in grid, z-direction */
h5_int64_t num_gridpx_m; /**< Read in number of points after 0(not counted here) in grid, r-direction*/
h5_int64_t num_gridpy_m; /**< Read in number of points after 0(not counted here) in grid, r-direction*/
h5_int64_t num_gridpz_m; /**< Read in number of points after 0(not counted here) in grid, z-direction*/
std::vector<double> FieldstrengthHz_m; /**< 3D array with Hz */
std::vector<double> FieldstrengthHx_m; /**< 3D array with Hx */
std::vector<double> FieldstrengthHy_m; /**< 3D array with Hy */
bool swap_m;
friend class Fieldmap;
friend class FM3DH5BlockBase;
};
inline bool FM3DH5Block_nonscale::isInside(const Vector_t &r) const
{
return ((r(0) >= xbegin_m && r(0) < xend_m) &&
(r(1) >= ybegin_m && r(1) < yend_m) &&
(r(2) >= zbegin_m && r(2) < zend_m));
}
#endif
#endif
\ No newline at end of file
// vi: set et ts=4 sw=4 sts=4:
// Local Variables:
// mode:c
// c-basic-offset: 4
// indent-tabs-mode: nil
// require-final-newline: nil
// End:
#ifndef CLASSIC_FIELDMAP3DMAGNETOSTATICH5BLOCK_HH
#define CLASSIC_FIELDMAP3DMAGNETOSTATICH5BLOCK_HH
#include "Fields/Fieldmap.h"
#include "hdf5.h"
#include "H5hut.h"
class FM3DMagnetoStaticH5Block: public Fieldmap {
//
// Class for magneto-static 3D field-maps stored in H5hut files.
//
// Copyright (c) 2020, Achim Gsell, Paul Scherrer Institut, Villigen PSI, Switzerland
// All rights reserved.
//
// This file is part of OPAL.
//
// OPAL is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// You should have received a copy of the GNU General Public License
// along with OPAL. If not, see <https://www.gnu.org/licenses/>.
//
#ifndef CLASSIC_FIELDMAP3DMAGNETOSTATICH5BLOCK_H
#define CLASSIC_FIELDMAP3DMAGNETOSTATICH5BLOCK_H
#include "Fields/FM3DH5BlockBase.h"
#include <vector>
class FM3DMagnetoStaticH5Block: public FM3DH5BlockBase {
public:
virtual bool getFieldstrength(const Vector_t &R, Vector_t &E, Vector_t &B) const;
virtual void getFieldDimensions(double &zBegin, double &zEnd, double &rBegin, double &rEnd) const;
virtual void getFieldDimensions(double &xIni, double &xFinal, double &yIni, double &yFinal, double &zIni, double &zFinal) const;
virtual bool getFieldDerivative(const Vector_t &R, Vector_t &E, Vector_t &B, const DiffDirection &dir) const;
virtual void swap();
virtual void getInfo(Inform *msg);
virtual double getFrequency() const;
virtual void setFrequency(double freq);
virtual bool isInside(const Vector_t &r) const;
virtual bool getFieldstrength (
const Vector_t &R, Vector_t &E, Vector_t &B) const;
private:
FM3DMagnetoStaticH5Block(std::string aFilename);
~FM3DMagnetoStaticH5Block();
virtual void readMap();
virtual void freeMap();
std::vector<h5_float64_t> FieldstrengthEz_m; /**< 3D array with Ez */
std::vector<h5_float64_t> FieldstrengthEx_m; /**< 3D array with Ex */
std::vector<h5_float64_t> FieldstrengthEy_m; /**< 3D array with Ey */
std::vector<h5_float64_t> FieldstrengthBz_m; /**< 3D array with Bz */
std::vector<h5_float64_t> FieldstrengthBx_m; /**< 3D array with Bx */
std::vector<h5_float64_t> FieldstrengthBy_m; /**< 3D array with By */
FM3DMagnetoStaticH5Block (
std::string aFilename);
h5_float64_t frequency_m;
virtual ~FM3DMagnetoStaticH5Block (
);
h5_float64_t xbegin_m;
h5_float64_t xend_m;
// int xcentral_idx_m;
virtual void readMap (
);
h5_float64_t ybegin_m;
h5_float64_t yend_m;
// int ycentral_idx_m;
virtual void freeMap (
);
h5_float64_t zbegin_m;
h5_float64_t zend_m;
virtual double getFrequency (
) const;
h5_float64_t hx_m; /**< length between points in grid, x-direction */
h5_float64_t hy_m; /**< length between points in grid, y-direction */
h5_float64_t hz_m; /**< length between points in grid, z-direction */
h5_int64_t num_gridpx_m; /**< Read in number of points after 0(not counted here) in grid, r-direction*/
h5_int64_t num_gridpy_m; /**< Read in number of points after 0(not counted here) in grid, r-direction*/
h5_int64_t num_gridpz_m; /**< Read in number of points after 0(not counted here) in grid, z-direction*/
std::vector<double> FieldstrengthBz_m; /**< 3D array with Hz */
std::vector<double> FieldstrengthBx_m; /**< 3D array with Hx */
std::vector<double> FieldstrengthBy_m; /**< 3D array with Hy */
bool swap_m;
friend class Fieldmap;
friend class FM3DH5BlockBase;
};
inline bool FM3DMagnetoStaticH5Block::isInside(const Vector_t &r) const
{
return ((r(0) >= xbegin_m && r(0) < xend_m) &&
(r(1) >= ybegin_m && r(1) < yend_m) &&
(r(2) >= zbegin_m && r(2) < zend_m));
}
#endif
#endif
\ No newline at end of file
// vi: set et ts=4 sw=4 sts=4:
// Local Variables:
// mode:c
// c-basic-offset: 4
// indent-tabs-mode: nil
// require-final-newline: nil
// End:
#include "Fields/Fieldmap.h"
#include "Fields/FM3DDynamic.h"
#include "Fields/FM3DH5BlockBase.h"
#include "Fields/FM3DH5Block.h"
#include "Fields/FM3DH5Block_nonscale.h"
#include "Fields/FM3DMagnetoStaticH5Block.h"
......
......@@ -103,6 +103,7 @@ public:
virtual void freeMap() = 0;
protected:
Fieldmap () {};
Fieldmap(const std::string& aFilename)
: Filename_m(aFilename),
lines_read_m(0),
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment