![]() |
![]() |
#include <MagneticField/Interpolation/src/TrapezoidalCartesianMFGrid.h>
Public Member Functions | |
void | dump () const |
virtual LocalPoint | fromGridFrame (double a, double b, double c) const |
find grid coordinates for point. For debugging and validation only. | |
virtual void | toGridFrame (const LocalPoint &p, double &a, double &b, double &c) const |
find grid coordinates for point. For debugging and validation only. | |
TrapezoidalCartesianMFGrid (binary_ifstream &istr, const GloballyPositioned< float > &vol) | |
virtual LocalVector | uncheckedValueInTesla (const LocalPoint &p) const |
Interpolated field value at given point; does not check for exceptions. | |
Private Attributes | |
bool | convertToLocal |
bool | increasingAlongX |
Trapezoid2RectangleMappingX | mapping_ |
The grid must have uniform spacing in two coordinates and increasing spacing in the other. Increasing spacing is supported only for x and y for the time being
Definition at line 20 of file TrapezoidalCartesianMFGrid.h.
TrapezoidalCartesianMFGrid::TrapezoidalCartesianMFGrid | ( | binary_ifstream & | istr, | |
const GloballyPositioned< float > & | vol | |||
) |
Definition at line 18 of file TrapezoidalCartesianMFGrid.cc.
References a, a1, b, convertToLocal, GenMuonPlsPt100GeV_cfg::cout, Vector3DBase< T, FrameTag >::dot(), dump(), lat::endl(), MFGrid::frame(), MFGrid3D::grid_, h, increasingAlongX, Grid1D< T >::lower(), mapping_, submitDQMOfflineCAF::nLines, offset, Trapezoid2RectangleMappingX::rectangle(), GloballyPositioned< T >::toGlobal(), GloballyPositioned< T >::toLocal(), Grid1D< T >::upper(), PV3DBase< T, VectorTag, FrameTag >::x(), PV3DBase< T, PVType, FrameType >::x(), x, PV3DBase< T, VectorTag, FrameTag >::y(), PV3DBase< T, PVType, FrameType >::y(), y, PV3DBase< T, PVType, FrameType >::z(), z, and PV3DBase< T, VectorTag, FrameTag >::z().
00020 : MFGrid3D(vol), increasingAlongX(false), convertToLocal(true) 00021 { 00022 00023 // The parameters read from the data files are given in global coordinates. 00024 // In version 85l, local frame has the same orientation of global frame for the reference 00025 // volume, i.e. the r.f. transformation is only a translation. 00026 // There is therefore no need to convert the field values to local coordinates. 00027 00028 // Check orientation of local reference frame: 00029 GlobalVector localXDir(frame().toGlobal(LocalVector(1,0,0))); 00030 GlobalVector localYDir(frame().toGlobal(LocalVector(0,1,0))); 00031 00032 if (localXDir.dot(GlobalVector(1,0,0)) > 0.999999 && 00033 localYDir.dot(GlobalVector(0,1,0)) > 0.999999) { 00034 // "null" rotation - requires no conversion... 00035 convertToLocal = false; 00036 } else if (localXDir.dot(GlobalVector(0,1,0)) > 0.999999 && 00037 localYDir.dot(GlobalVector(1,0,0)) > 0.999999) { 00038 // Typical orientation if master volume is in sector 1 00039 convertToLocal = true; 00040 } else { 00041 convertToLocal = true; 00042 // Nothing wrong in principle, but this is not expected 00043 cout << "WARNING: TrapezoidalCartesianMFGrid: unexpected orientation: x: " 00044 << localXDir << " y: " << localYDir << endl; 00045 } 00046 00047 int n1, n2, n3; 00048 inFile >> n1 >> n2 >> n3; 00049 double xref, yref, zref; 00050 inFile >> xref >> yref >> zref; 00051 double step1, step2, step3; 00052 inFile >> step1 >> step2 >> step3; 00053 00054 double BasicDistance1[3][3]; // linear step 00055 double BasicDistance2[3][3]; // linear offset 00056 bool easya, easyb, easyc; 00057 00058 inFile >> BasicDistance1[0][0] >> BasicDistance1[1][0] >> BasicDistance1[2][0]; 00059 inFile >> BasicDistance1[0][1] >> BasicDistance1[1][1] >> BasicDistance1[2][1]; 00060 inFile >> BasicDistance1[0][2] >> BasicDistance1[1][2] >> BasicDistance1[2][2]; 00061 inFile >> BasicDistance2[0][0] >> BasicDistance2[1][0] >> BasicDistance2[2][0]; 00062 inFile >> BasicDistance2[0][1] >> BasicDistance2[1][1] >> BasicDistance2[2][1]; 00063 inFile >> BasicDistance2[0][2] >> BasicDistance2[1][2] >> BasicDistance2[2][2]; 00064 inFile >> easya >> easyb >> easyc; 00065 00066 vector<BVector> fieldValues; 00067 float Bx, By, Bz; 00068 int nLines = n1*n2*n3; 00069 for (int iLine=0; iLine<nLines; ++iLine){ 00070 inFile >> Bx >> By >> Bz; 00071 if (convertToLocal) { 00072 // Preserve double precision! 00073 Vector3DBase<double, LocalTag> lB = frame().toLocal(Vector3DBase<double, GlobalTag>(Bx,By,Bz)); 00074 fieldValues.push_back(BVector(lB.x(), lB.y(), lB.z())); 00075 00076 } else { 00077 fieldValues.push_back(BVector(Bx,By,Bz)); 00078 } 00079 } 00080 // check completeness 00081 string lastEntry; 00082 inFile >> lastEntry; 00083 if (lastEntry != "complete") { 00084 cout << "ERROR during file reading: file is not complete" << endl; 00085 } 00086 00087 // In version 1103l the reference sector is at phi=0 so that local y is along global X. 00088 // The increasing grid steps for the same volume are then along Y instead than along X. 00089 // To use Trapezoid2RectangleMappingX to map the trapezoidal geometry to a rectangular 00090 // cartesian geometry, we have to exchange (global) X and Y appropriately when constructing 00091 // it. 00092 int nx, ny; 00093 double stepx, stepy, stepz; 00094 double dstep, offset; 00095 if (!easya && easyb && easyc) { 00096 // Increasing grid spacing is on x 00097 increasingAlongX = true; 00098 nx = n1; 00099 ny = n2; 00100 stepx = step1; 00101 stepy = step2; 00102 stepz = step3; 00103 dstep = BasicDistance1[0][1]; 00104 offset = BasicDistance2[0][1]; 00105 } else if (easya && !easyb && easyc) { 00106 // Increasing grid spacing is on y 00107 increasingAlongX = false; 00108 nx = n2; 00109 ny = n1; 00110 stepx = step2; 00111 stepy = step1; 00112 stepz = -step3; 00113 dstep = BasicDistance1[1][0]; 00114 offset = BasicDistance2[1][0]; 00115 00116 } else { 00117 // Increasing spacing on z or on > 1 coordinate not supported 00118 throw MagGeometryError("TrapezoidalCartesianMFGrid only implemented for first or second coordinate"); 00119 } 00120 00121 double a = stepx * (nx -1); // first base 00122 double b = a + dstep * (ny-1) * (nx-1); // second base 00123 double h = stepy * (ny-1); // height 00124 double delta = -offset * (ny-1); // offset between two bases 00125 double baMinus1 = dstep*(ny-1) / stepx; // (b*a) - 1 00126 00127 GlobalPoint grefp( xref, yref, zref); 00128 LocalPoint lrefp = frame().toLocal( grefp); 00129 00130 if (fabs(baMinus1) > 0.000001) { 00131 double b_over_a = 1 + baMinus1; 00132 double a1 = delta/baMinus1; 00133 00134 #ifdef DEBUG_GRID 00135 cout << "Trapeze size (a,b,h) = " << a << "," << b << "," << h << endl; 00136 cout << "Global origin " << grefp << endl; 00137 cout << "Local origin " << lrefp << endl; 00138 cout << "a1 = " << a1 << endl; 00139 #endif 00140 00141 // FIXME ASSUMPTION: here we assume that the local reference frame is oriented with X along 00142 // the direction of where the grid is not uniform. This is the case for all current geometries 00143 double x0 = lrefp.x() + a1; 00144 double y0 = lrefp.y() + h/2.; 00145 mapping_ = Trapezoid2RectangleMappingX( x0, y0, b_over_a, h); 00146 } 00147 else { // parallelogram 00148 mapping_ = Trapezoid2RectangleMappingX( 0, 0, delta/h); 00149 } 00150 00151 // transform reference point to grid frame 00152 double xrec, yrec; 00153 mapping_.rectangle( lrefp.x(), lrefp.y(), xrec, yrec); 00154 00155 Grid1D<double> gridX( xrec, xrec + (a+b)/2., nx); 00156 Grid1D<double> gridY( yrec, yrec + h, ny); 00157 Grid1D<double> gridZ( lrefp.z(), lrefp.z() + stepz*(n3-1), n3); 00158 00159 #ifdef DEBUG_GRID 00160 cout << " GRID X range: local " << gridX.lower() << " - " << gridX.upper() 00161 <<" global: " << (frame().toGlobal(LocalPoint(gridX.lower(),0,0))).y() << " - " 00162 << (frame().toGlobal(LocalPoint(gridX.upper(),0,0))).y() << endl; 00163 00164 cout << " GRID Y range: local " << gridY.lower() << " - " << gridY.upper() 00165 << " global: " << (frame().toGlobal(LocalPoint(0,gridY.lower(),0))).x() << " - " 00166 << (frame().toGlobal(LocalPoint(0,gridY.upper(),0))).x() << endl; 00167 00168 cout << " GRID Z range: local " << gridZ.lower() << " - " << gridZ.upper() 00169 << " global: " << (frame().toGlobal(LocalPoint(0,0,gridZ.lower()))).z() << " " 00170 << (frame().toGlobal(LocalPoint(0,0,gridZ.upper()))).z() << endl; 00171 #endif 00172 00173 if (increasingAlongX) { 00174 grid_ = GridType( gridX, gridY, gridZ, fieldValues); 00175 } else { 00176 // The reason why gridY and gridX have to be exchanged is because Grid3D::index(i,j,k) 00177 // assumes a specific order for the fieldValues, and we cannot rearrange this vector. 00178 // Given that we exchange grids, we will have to exchange the outpouts of mapping_rectangle() 00179 // and the inputs of mapping_.trapezoid() in the following... 00180 grid_ = GridType( gridY, gridX, gridZ, fieldValues); 00181 } 00182 00183 00184 #ifdef DEBUG_GRID 00185 dump(); 00186 #endif 00187 }
Reimplemented from MFGrid.
Definition at line 189 of file TrapezoidalCartesianMFGrid.cc.
References GenMuonPlsPt100GeV_cfg::cout, Grid3D< Value, T >::data(), lat::endl(), MFGrid3D::grid_, Grid3D< Value, T >::grida(), Grid3D< Value, T >::gridb(), and Grid3D< Value, T >::gridc().
Referenced by TrapezoidalCartesianMFGrid().
00190 { 00191 00192 cout << endl << "Dump of TrapezoidalCartesianMFGrid" << endl; 00193 // cout << "Number of points from file " 00194 // << n1 << " " << n2 << " " << n3 << endl; 00195 cout << "Number of points from Grid1D " 00196 << grid_.grida().nodes() << " " << grid_.gridb().nodes() << " " << grid_.gridc().nodes() << endl; 00197 00198 // cout << "Reference Point from file " 00199 // << xref << " " << yref << " " << zref << endl; 00200 cout << "Reference Point from Grid1D " 00201 << grid_.grida().lower() << " " << grid_.gridb().lower() << " " << grid_.gridc().lower() << endl; 00202 00203 // cout << "Basic Distance from file " 00204 // << stepx << " " << stepy << " " << stepz << endl; 00205 cout << "Basic Distance from Grid1D " 00206 << grid_.grida().step() << " " << grid_.gridb().step() << " " << grid_.gridc().step() << endl; 00207 00208 00209 cout << "Dumping " << grid_.data().size() << " field values " << endl; 00210 // grid_.dump(); 00211 00212 00213 // Dump ALL grid points and values 00214 // CAVEAT: if convertToLocal = true in the ctor, points have been converted to LOCAL 00215 // coordinates. To match those from .table files they have to be converted back to global 00216 // for (int j=0; j < grid_.gridb().nodes(); j++) { 00217 // for (int i=0; i < grid_.grida().nodes(); i++) { 00218 // for (int k=0; k < grid_.gridc().nodes(); k++) { 00219 // cout << i << " " << j << " " << k << " " 00220 // << frame().toGlobal(LocalPoint(nodePosition(i,j,k))) << " " 00221 // << nodeValue(i,j,k) << endl; 00222 // } 00223 // } 00224 // } 00225 00226 00227 }
MFGrid::LocalPoint TrapezoidalCartesianMFGrid::fromGridFrame | ( | double | a, | |
double | b, | |||
double | c | |||
) | const [virtual] |
find grid coordinates for point. For debugging and validation only.
Implements MFGrid.
Definition at line 264 of file TrapezoidalCartesianMFGrid.cc.
References increasingAlongX, mapping_, and Trapezoid2RectangleMappingX::trapezoid().
00265 { 00266 double xtrap, ytrap; 00267 if (increasingAlongX) { 00268 mapping_.trapezoid( a, b, xtrap, ytrap); 00269 } else { 00270 mapping_.trapezoid( b, a, xtrap, ytrap); 00271 } 00272 return LocalPoint( xtrap, ytrap, c); 00273 }
void TrapezoidalCartesianMFGrid::toGridFrame | ( | const LocalPoint & | p, | |
double & | a, | |||
double & | b, | |||
double & | c | |||
) | const [virtual] |
find grid coordinates for point. For debugging and validation only.
Implements MFGrid.
Definition at line 252 of file TrapezoidalCartesianMFGrid.cc.
References increasingAlongX, mapping_, Trapezoid2RectangleMappingX::rectangle(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().
00254 { 00255 if (increasingAlongX) { 00256 mapping_.rectangle( p.x(), p.y(), a, b); 00257 } else { 00258 mapping_.rectangle( p.x(), p.y(), b, a); 00259 } 00260 00261 c = p.z(); 00262 }
MFGrid::LocalVector TrapezoidalCartesianMFGrid::uncheckedValueInTesla | ( | const LocalPoint & | p | ) | const [virtual] |
Interpolated field value at given point; does not check for exceptions.
Implements MFGrid3D.
Definition at line 229 of file TrapezoidalCartesianMFGrid.cc.
References MFGrid3D::grid_, increasingAlongX, LinearGridInterpolator3D< Value, T >::interpolate(), mapping_, Trapezoid2RectangleMappingX::rectangle(), std::swap(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().
00230 { 00231 00232 double xrec, yrec; 00233 mapping_.rectangle( p.x(), p.y(), xrec, yrec); 00234 00235 00236 // cout << "TrapezoidalCartesianMFGrid::valueInTesla at local point " << p << endl; 00237 // cout << p.x() << " " << p.y() 00238 // << " transformed to grid frame: " << xrec << " " << yrec << endl; 00239 00240 LinearGridInterpolator3D<GridType::ValueType, GridType::Scalar> interpol( grid_); 00241 00242 if (!increasingAlongX) { 00243 std::swap(xrec,yrec); 00244 // B values are already converted to local coord!!! otherwise we should: 00245 // GridType::ValueType value = interpol.interpolate( xrec, yrec, p.z()); 00246 // return LocalVector(value.y(),value.x(),-value.z()); 00247 } 00248 00249 return LocalVector(interpol.interpolate( xrec, yrec, p.z())); 00250 }
Definition at line 38 of file TrapezoidalCartesianMFGrid.h.
Referenced by TrapezoidalCartesianMFGrid().
Definition at line 37 of file TrapezoidalCartesianMFGrid.h.
Referenced by fromGridFrame(), toGridFrame(), TrapezoidalCartesianMFGrid(), and uncheckedValueInTesla().
Definition at line 36 of file TrapezoidalCartesianMFGrid.h.
Referenced by fromGridFrame(), toGridFrame(), TrapezoidalCartesianMFGrid(), and uncheckedValueInTesla().