diff --git a/CMakeLists.txt b/CMakeLists.txt index 96e0c24753939..0abf370d785ba 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -163,6 +163,7 @@ add_subdirectory (Data) add_subdirectory (Generators) add_subdirectory (its) add_subdirectory (passive) +add_subdirectory (MathUtils) add_subdirectory (field) add_subdirectory (devices) Add_Subdirectory(macro) diff --git a/MathUtils/CMakeLists.txt b/MathUtils/CMakeLists.txt new file mode 100644 index 0000000000000..c636b7c4bead3 --- /dev/null +++ b/MathUtils/CMakeLists.txt @@ -0,0 +1,32 @@ +# Create a library called "MathUtils" which includes the source files given in +# the array. +# The extension is already found. Any number of sources could be listed here. + +set(INCLUDE_DIRECTORIES +${CMAKE_SOURCE_DIR}/MathUtils +${BASE_INCLUDE_DIRECTORIES} +${ROOT_INCLUDE_DIR} +) + +include_directories( ${INCLUDE_DIRECTORIES}) + +set(LINK_DIRECTORIES +${CMAKE_SOURCE_DIR}/MathUtils +${FAIRROOT_LIBRARY_DIR} +${ROOT_LIBRARY_DIR} +) + +link_directories( ${LINK_DIRECTORIES}) + +set(SRCS +Chebyshev3D.cxx +Chebyshev3DCalc.cxx +) + +Set(HEADERS) +Set(LINKDEF MathUtilsLinkDef.h) +Set(LIBRARY_NAME MathUtils) +Set(DEPENDENCIES Cint Core) + +GENERATE_LIBRARY() + diff --git a/field/Chebyshev3D.cxx b/MathUtils/Chebyshev3D.cxx similarity index 93% rename from field/Chebyshev3D.cxx rename to MathUtils/Chebyshev3D.cxx index 4881c8eba63ab..1e34a588f7d00 100644 --- a/field/Chebyshev3D.cxx +++ b/MathUtils/Chebyshev3D.cxx @@ -14,13 +14,15 @@ #include "Chebyshev3DCalc.h" #include "FairLogger.h" -using namespace AliceO2::Field; +using namespace AliceO2::MathUtils; ClassImp(Chebyshev3D) +const Float_t Chebyshev3D::sMinimumPrecision = 1.e-12f; + Chebyshev3D::Chebyshev3D() : mOutputArrayDimension(0), - mPrecision(0), + mPrecision(sMinimumPrecision), mChebyshevParameter(1), mMaxCoefficients(0), mTemporaryUserResults(0), @@ -111,11 +113,11 @@ Chebyshev3D::Chebyshev3D(FILE* stream) } #ifdef _INC_CREATION_Chebyshev3D_ -Chebyshev3D::Chebyshev3D(const char* funName, int DimOut, const Float_t* bmin, const Float_t* bmax, Int_t* npoints, - Float_t prec) +Chebyshev3D::Chebyshev3D(const char* funName, int dimOut, const Float_t* bmin, const Float_t* bmax, + const Int_t* npoints, Float_t prec, const Float_t* precD) : TNamed(funName, funName), mOutputArrayDimension(0), - mPrecision(TMath::Max(1.E-12f, prec)), + mPrecision(TMath::Max(sMinimumPrecision, prec)), mChebyshevParameter(1), mMaxCoefficients(0), mTemporaryUserResults(0), @@ -124,7 +126,7 @@ Chebyshev3D::Chebyshev3D(const char* funName, int DimOut, const Float_t* bmin, c mUserMacro(0), mLogger(FairLogger::GetLogger()) { - if (DimOut < 1) { + if (dimOut < 1) { Error("Chebyshev3D", "Requested output dimension is %d\nStop\n", mOutputArrayDimension); exit(1); } @@ -134,7 +136,7 @@ Chebyshev3D::Chebyshev3D(const char* funName, int DimOut, const Float_t* bmin, c mTemporaryChebyshevGridOffs[i] = 0.; mTemporaryCoefficient[i] = 0; } - setDimOut(DimOut); + setDimOut(dimOut,precD); prepareBoundaries(bmin, bmax); defineGrid(npoints); setuserFunction(funName); @@ -143,10 +145,10 @@ Chebyshev3D::Chebyshev3D(const char* funName, int DimOut, const Float_t* bmin, c #endif #ifdef _INC_CREATION_Chebyshev3D_ -Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, Float_t* bmax, Int_t* npoints, - Float_t prec) +Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int dimOut, const Float_t* bmin, const Float_t* bmax, + const Int_t* npoints, Float_t prec, const Float_t* precD) : mOutputArrayDimension(0), - mPrecision(TMath::Max(1.E-12f, prec)), + mPrecision(TMath::Max(sMinimumPrecision, prec)), mChebyshevParameter(1), mMaxCoefficients(0), mTemporaryUserResults(0), @@ -155,11 +157,7 @@ Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, mUserMacro(0), mLogger(FairLogger::GetLogger()) { - if (DimOut < 1) { - Error("Chebyshev3D", "Requested output dimension is %d\nStop\n", mOutputArrayDimension); - exit(1); - } - if (DimOut < 1) { + if (dimOut < 1) { Error("Chebyshev3D", "Requested output dimension is %d\nStop\n", mOutputArrayDimension); exit(1); } @@ -169,7 +167,7 @@ Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, mTemporaryChebyshevGridOffs[i] = 0.; mTemporaryCoefficient[i] = 0; } - setDimOut(DimOut); + setDimOut(dimOut,precD); prepareBoundaries(bmin, bmax); defineGrid(npoints); setuserFunction(ptr); @@ -178,10 +176,10 @@ Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, #endif #ifdef _INC_CREATION_Chebyshev3D_ -Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, Float_t* bmax, Int_t* npX, Int_t* npY, - Int_t* npZ, Float_t prec) +Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int dimOut, const Float_t* bmin, const Float_t* bmax, + const Int_t* npX, const Int_t* npY, const Int_t* npZ, Float_t prec, const Float_t* precD) : mOutputArrayDimension(0), - mPrecision(TMath::Max(1.E-12f, prec)), + mPrecision(TMath::Max(sMinimumPrecision, prec)), mChebyshevParameter(1), mMaxCoefficients(0), mTemporaryUserResults(0), @@ -190,11 +188,7 @@ Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, mUserMacro(0), mLogger(FairLogger::GetLogger()) { - if (DimOut < 1) { - Error("Chebyshev3D", "Requested output dimension is %d\nStop\n", mOutputArrayDimension); - exit(1); - } - if (DimOut < 1) { + if (dimOut < 1) { Error("Chebyshev3D", "Requested output dimension is %d\nStop\n", mOutputArrayDimension); exit(1); } @@ -204,7 +198,7 @@ Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, mTemporaryChebyshevGridOffs[i] = 0.; mTemporaryCoefficient[i] = 0; } - setDimOut(DimOut); + setDimOut(dimOut,precD); prepareBoundaries(bmin, bmax); setuserFunction(ptr); @@ -218,10 +212,10 @@ Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, #endif #ifdef _INC_CREATION_Chebyshev3D_ -Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, Float_t* bmax, Float_t prec, - Bool_t run) +Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int dimOut, const Float_t* bmin, const Float_t* bmax, + Float_t prec, Bool_t run, const Float_t* precD) : mOutputArrayDimension(0), - mPrecision(TMath::Max(1.E-12f, prec)), + mPrecision(TMath::Max(sMinimumPrecision, prec)), mChebyshevParameter(1), mMaxCoefficients(0), mTemporaryUserResults(0), @@ -230,11 +224,11 @@ Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, mUserMacro(0), mLogger(FairLogger::GetLogger()) { - if (DimOut != 3) { + if (dimOut != 3) { Error("Chebyshev3D", "This constructor works only for 3D fits, %dD fit was requested\n", mOutputArrayDimension); exit(1); } - if (DimOut < 1) { + if (dimOut < 1) { Error("Chebyshev3D", "Requested output dimension is %d\nStop\n", mOutputArrayDimension); exit(1); } @@ -244,7 +238,7 @@ Chebyshev3D::Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, mTemporaryChebyshevGridOffs[i] = 0.; mTemporaryCoefficient[i] = 0; } - setDimOut(DimOut); + setDimOut(dimOut, precD); prepareBoundaries(bmin, bmax); setuserFunction(ptr); @@ -458,7 +452,7 @@ Int_t Chebyshev3D::calculateChebyshevCoefficients(const Float_t* funval, int np, #endif #ifdef _INC_CREATION_Chebyshev3D_ -void Chebyshev3D::defineGrid(Int_t* npoints) +void Chebyshev3D::defineGrid(const Int_t* npoints) { // prepare the grid of Chebyshev roots in each dimension const int kMinPoints = 1; @@ -519,8 +513,6 @@ Int_t Chebyshev3D::chebyshevFit(int dmOut) Float_t* tmpCoef2D = new Float_t[mNumberOfPoints[0] * mNumberOfPoints[1]]; Float_t* tmpCoef1D = new Float_t[maxDim]; - Float_t rTiny = 0.1 * mPrecision / Float_t(maxDim); // neglect coefficient below this threshold - // 1D Cheb.fit for 0-th dimension at current steps of remaining dimensions int ncmax = 0; @@ -528,6 +520,10 @@ Int_t Chebyshev3D::chebyshevFit(int dmOut) fflush(stdout); Chebyshev3DCalc* cheb = getChebyshevCalc(dmOut); + Float_t prec = cheb->getPrecision(); + if (precinitRows(nRows); // create needed arrays; + cheb->initializeRows(nRows); // create needed arrays; UShort_t* nColsAtRow = cheb->getNumberOfColumnsAtRow(); - UShort_t* colAtRowBg = cheb->GetColAtRowBg(); + UShort_t* colAtRowBg = cheb->getColAtRowBg(); int nCols = 0; int nElemBound2D = 0; for (int id0 = 0; id0 < nRows; id0++) { @@ -651,12 +647,12 @@ Int_t Chebyshev3D::chebyshevFit(int dmOut) nCols = nColsAtRow[id0]; } } - cheb->initCols(nCols); + cheb->initializeColumns(nCols); delete[] tmpCols; // create the 2D matrix defining the boundary of significance for 3D coeffs.matrix // and count the number of siginifacnt coefficients - cheb->InitElemBound2D(nElemBound2D); + cheb->initializeElementBound2D(nElemBound2D); UShort_t* coefBound2D0 = cheb->getCoefficientBound2D0(); UShort_t* coefBound2D1 = cheb->getCoefficientBound2D1(); mMaxCoefficients = 0; // redefine number of coeffs @@ -821,7 +817,7 @@ void Chebyshev3D::loadData(FILE* stream) } } -void Chebyshev3D::setDimOut(const int d) +void Chebyshev3D::setDimOut(const int d, const float* prec) { // init output dimensions mOutputArrayDimension = d; @@ -831,6 +827,8 @@ void Chebyshev3D::setDimOut(const int d) mTemporaryUserResults = new Float_t[mOutputArrayDimension]; mChebyshevParameter.Delete(); for (int i = 0; i < d; i++) { + Chebyshev3DCalc* clc = new Chebyshev3DCalc(); + clc->setPrecision(prec && prec[i]>sMinimumPrecision ? prec[i] : mPrecision); mChebyshevParameter.AddAtAndExpand(new Chebyshev3DCalc(), i); } } @@ -860,6 +858,10 @@ TH1* Chebyshev3D::TestRMS(int idim, int npoints, TH1* histo) if (!histo) { histo = new TH1D(GetName(), "Control: Function - Parametrization", 100, -2 * mPrecision, 2 * mPrecision); } + + float prc = getChebyshevCalc(idim)->getPrecision(); + if (prcRndmArray(3, (Float_t*)mTemporaryCoefficient); for (int i = 3; i--;) { @@ -877,7 +879,7 @@ TH1* Chebyshev3D::TestRMS(int idim, int npoints, TH1* histo) #ifdef _INC_CREATION_Chebyshev3D_ -void Chebyshev3D::estimateNumberOfPoints(float Prec, int gridBC[3][3], Int_t npd1, Int_t npd2, Int_t npd3) +void Chebyshev3D::estimateNumberOfPoints(float prec, int gridBC[3][3], Int_t npd1, Int_t npd2, Int_t npd3) { // Estimate number of points to generate a training data const int kScp = 9; @@ -904,7 +906,7 @@ void Chebyshev3D::estimateNumberOfPoints(float Prec, int gridBC[3][3], Int_t npd xyz[id1] = mMinBoundaries[id1] + kScl[i1] * (mMaxBoundaries[id1] - mMinBoundaries[id1]); for (int i2 = 0; i2 < kScp; i2++) { xyz[id2] = mMinBoundaries[id2] + kScl[i2] * (mMaxBoundaries[id2] - mMinBoundaries[id2]); - int* npt = getNcNeeded(xyz, idim, dimMN, dimMX, Prec, npdTst[idim]); // npoints for Bx,By,Bz + int* npt = getNcNeeded(xyz, idim, dimMN, dimMX, prec, npdTst[idim]); // npoints for Bx,By,Bz for (int ib = 0; ib < 3; ib++) { if (npt[ib] > gridBC[ib][idim]) { gridBC[ib][idim] = npt[ib]; diff --git a/field/Chebyshev3D.h b/MathUtils/Chebyshev3D.h similarity index 90% rename from field/Chebyshev3D.h rename to MathUtils/Chebyshev3D.h index baa899d72e271..b12a578c3fa08 100644 --- a/field/Chebyshev3D.h +++ b/MathUtils/Chebyshev3D.h @@ -2,8 +2,8 @@ /// \brief Definition of the Cheb3D class /// \author ruben.shahoyan@cern.ch 09/09/2006 -#ifndef ALICEO2_FIELD_CHEBYSHEV3D_H_ -#define ALICEO2_FIELD_CHEBYSHEV3D_H_ +#ifndef ALICEO2_MATHUTILS_CHEBYSHEV3D_H_ +#define ALICEO2_MATHUTILS_CHEBYSHEV3D_H_ #include #include @@ -21,7 +21,7 @@ class FairLogger; namespace AliceO2 { -namespace Field { +namespace MathUtils { /// Chebyshev3D produces the interpolation of the user 3D->NDimOut arbitrary function supplied in /// "void (*fcn)(float* inp,float* out)" format either in a separate macro file or as a function pointer. Only @@ -69,9 +69,10 @@ class Chebyshev3D : public TNamed { /// \param bmax : array of 3 elements with the upper boundaries of the region where the function is defined /// \param npoints : array of 3 elements with the number of points to compute in each of 3 dimension /// \param prec : max allowed absolute difference between the user function and computed parameterization on the - /// requested grid - Chebyshev3D(const char* funName, Int_t DimOut, const Float_t* bmin, const Float_t* bmax, Int_t* npoints, - Float_t prec = 1E-6); + /// requested grid, common for all 1D components + /// \param precD : optional precison per component + Chebyshev3D(const char* funName, Int_t dimOut, const Float_t* bmin, const Float_t* bmax, const Int_t* npoints, + Float_t prec = 1E-6, const Float_t* precD=0); /// Construct the parameterization for the function /// \param ptr : pointer on the function: void fun(Float_t * inp,Float_t * out) /// \param DimOut : dimension of the vector computed by the user function @@ -79,9 +80,10 @@ class Chebyshev3D : public TNamed { /// \param bmax : array of 3 elements with the upper boundaries of the region where the function is defined /// \param npoints : array of 3 elements with the number of points to compute in each of 3 dimension /// \param prec : max allowed absolute difference between the user function and computed parameterization on the - /// requested grid - Chebyshev3D(void (*ptr)(float*, float*), Int_t DimOut, Float_t* bmin, Float_t* bmax, Int_t* npoints, - Float_t prec = 1E-6); + /// requested grid, common for all 1D components + /// \param precD : optional precison per component + Chebyshev3D(void (*ptr)(float*, float*), Int_t dimOut, const Float_t* bmin, const Float_t* bmax, const Int_t* npoints, + Float_t prec = 1E-6, const Float_t* precD=0); /// Construct very economic parameterization for the function /// \param ptr : pointer on the function: void fun(Float_t * inp,Float_t * out) /// \param DimOut : dimension of the vector computed by the user function @@ -91,18 +93,21 @@ class Chebyshev3D : public TNamed { /// \param npY : array of 3 elements with the number of points to compute in each dimension for 2nd component /// \param npZ : array of 3 elements with the number of points to compute in each dimension for 3d component /// \param prec : max allowed absolute difference between the user function and computed parameterization on the - /// requested grid - Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, Float_t* bmax, Int_t* npX, Int_t* npY, Int_t* npZ, - Float_t prec = 1E-6); + /// requested grid, common for all 1D components + /// \param precD : optional precison per component + Chebyshev3D(void (*ptr)(float*, float*), int dimOut, const Float_t* bmin, const Float_t* bmax, + const Int_t* npX, const Int_t* npY, const Int_t* npZ, + Float_t prec = 1E-6, const Float_t* precD=0); /// Construct very economic parameterization for the function with automatic calculation of the root's grid /// \param ptr : pointer on the function: void fun(Float_t * inp,Float_t * out) /// \param DimOut : dimension of the vector computed by the user function /// \param bmin : array of 3 elements with the lower boundaries of the region where the function is defined /// \param bmax : array of 3 elements with the upper boundaries of the region where the function is defined /// \param prec : max allowed absolute difference between the user function and computed parameterization on the - /// \param requested grid - Chebyshev3D(void (*ptr)(float*, float*), int DimOut, Float_t* bmin, Float_t* bmax, Float_t prec = 1E-6, - Bool_t run = kTRUE); + /// \param requested grid, common for all 1D components + /// \param precD : optional precison per component + Chebyshev3D(void (*ptr)(float*, float*), int DimOut, const Float_t* bmin, const Float_t* bmax, Float_t prec = 1E-6, + Bool_t run = kTRUE, const Float_t* precD=0); #endif ~Chebyshev3D() @@ -163,8 +168,8 @@ class Chebyshev3D : public TNamed { #ifdef _INC_CREATION_Chebyshev3D_ void invertSign(); - int* getNcNeeded(float xyz[3], int DimVar, float mn, float mx, float prec, Int_t npCheck = 30); - void estimateNumberOfPoints(float Prec, int gridBC[3][3], Int_t npd1 = 30, Int_t npd2 = 30, Int_t npd3 = 30); + int* getNcNeeded(float xyz[3], int dimVar, float mn, float mx, float prec, Int_t npCheck = 30); + void estimateNumberOfPoints(float prec, int gridBC[3][3], Int_t npd1 = 30, Int_t npd2 = 30, Int_t npd3 = 30); void saveData(const char* outfile, Bool_t append = kFALSE) const; void saveData(FILE* stream = stdout) const; @@ -177,12 +182,12 @@ class Chebyshev3D : public TNamed { protected: void Clear(const Option_t* option = ""); - void setDimOut(const int d); + void setDimOut(const int d, const float* prec=0); void prepareBoundaries(const Float_t* bmin, const Float_t* bmax); #ifdef _INC_CREATION_Chebyshev3D_ void evaluateUserFunction(); - void defineGrid(Int_t* npoints); + void defineGrid(const Int_t* npoints); Int_t chebyshevFit(); // fit all output dimensions Int_t chebyshevFit(int dmOut); void setPrecision(float prec) @@ -221,7 +226,9 @@ class Chebyshev3D : public TNamed { TMethodCall* mUserMacro; //! Pointer to MethodCall for function from user macro FairLogger* mLogger; //! - ClassDef(AliceO2::Field::Chebyshev3D, 2) // Chebyshev parametrization for 3D->N function + static const Float_t sMinimumPrecision; ///< minimum precision allowed + + ClassDef(AliceO2::MathUtils::Chebyshev3D, 2) // Chebyshev parametrization for 3D->N function }; /// Checks if the point is inside of the fitted box diff --git a/field/Chebyshev3DCalc.cxx b/MathUtils/Chebyshev3DCalc.cxx similarity index 89% rename from field/Chebyshev3DCalc.cxx rename to MathUtils/Chebyshev3DCalc.cxx index 651b38caa7a9e..c201427d09286 100644 --- a/field/Chebyshev3DCalc.cxx +++ b/MathUtils/Chebyshev3DCalc.cxx @@ -6,7 +6,7 @@ #include #include "Chebyshev3DCalc.h" -using namespace AliceO2::Field; +using namespace AliceO2::MathUtils; ClassImp(Chebyshev3DCalc) @@ -15,6 +15,7 @@ Chebyshev3DCalc::Chebyshev3DCalc() mNumberOfRows(0), mNumberOfColumns(0), mNumberOfElementsBound2D(0), + mPrecision(0), mNumberOfColumnsAtRow(0), mColumnAtRowBeginning(0), mCoefficientBound2D0(0), @@ -31,6 +32,7 @@ Chebyshev3DCalc::Chebyshev3DCalc(const Chebyshev3DCalc& src) mNumberOfRows(src.mNumberOfRows), mNumberOfColumns(src.mNumberOfColumns), mNumberOfElementsBound2D(src.mNumberOfElementsBound2D), + mPrecision(src.mPrecision), mNumberOfColumnsAtRow(0), mColumnAtRowBeginning(0), mCoefficientBound2D0(0), @@ -82,6 +84,7 @@ Chebyshev3DCalc::Chebyshev3DCalc(FILE* stream) mNumberOfRows(0), mNumberOfColumns(0), mNumberOfElementsBound2D(0), + mPrecision(0), mNumberOfColumnsAtRow(0), mColumnAtRowBeginning(0), mCoefficientBound2D0(0), @@ -102,6 +105,7 @@ Chebyshev3DCalc& Chebyshev3DCalc::operator=(const Chebyshev3DCalc& rhs) mNumberOfCoefficients = rhs.mNumberOfCoefficients; mNumberOfRows = rhs.mNumberOfRows; mNumberOfColumns = rhs.mNumberOfColumns; + mPrecision = rhs.mPrecision; if (rhs.mNumberOfColumnsAtRow) { mNumberOfColumnsAtRow = new UShort_t[mNumberOfRows]; for (int i = mNumberOfRows; i--;) { @@ -176,7 +180,8 @@ void Chebyshev3DCalc::Clear(const Option_t*) void Chebyshev3DCalc::Print(const Option_t*) const { - printf("Chebyshev parameterization data %s for 3D->1 function.\n", GetName()); + printf("Chebyshev parameterization data %s for 3D->1 function, precision: %e\n", + GetName(),mPrecision); int nmax3d = 0; for (int i = mNumberOfElementsBound2D; i--;) { if (mCoefficientBound2D0[i] > nmax3d) { @@ -258,7 +263,7 @@ Float_t Chebyshev3DCalc::evaluateDerivative2(int dim1, int dim2, const Float_t* : chebyshevEvaluation1D(par[0], mTemporaryCoefficients1D, mNumberOfRows); } -#ifdef _INC_CREATION_ALICHEB3D_ +#ifdef _INC_CREATION_Chebyshev3D_ void Chebyshev3DCalc::saveData(const char* outfile, Bool_t append) const { TString strf = outfile; @@ -269,7 +274,7 @@ void Chebyshev3DCalc::saveData(const char* outfile, Bool_t append) const } #endif -#ifdef _INC_CREATION_ALICHEB3D_ +#ifdef _INC_CREATION_Chebyshev3D_ void Chebyshev3DCalc::saveData(FILE* stream) const { fprintf(stream, "#\nSTART %s\n", GetName()); @@ -289,6 +294,9 @@ void Chebyshev3DCalc::saveData(FILE* stream) const for (int i = 0; i < mNumberOfCoefficients; i++) { fprintf(stream, "%+.8e\n", mCoefficients[i]); } + fprintf(stream,"# Precision\n"); + fprintf(stream,"%+.8e\n",mPrecision); + // fprintf(stream, "END %s\n", GetName()); } #endif @@ -315,7 +323,7 @@ void Chebyshev3DCalc::loadData(FILE* stream) readLine(buffs, stream); // NRows mNumberOfRows = buffs.Atoi(); - if (mNumberOfRows < 1) { + if (mNumberOfRows < 0 || !buffs.IsDigit()) { Error("LoadData", "Expected: '', found \"%s\"\nStop\n", buffs.Data()); exit(1); } @@ -345,16 +353,15 @@ void Chebyshev3DCalc::loadData(FILE* stream) mNumberOfCoefficients += mCoefficientBound2D0[i]; } - if (mNumberOfCoefficients <= 0) { - Error("LoadData", "Negtive (%d) number of Chebychef coeffs. is obtained.\nStop\n", mNumberOfCoefficients); - exit(1); - } - initializeCoefficients(mNumberOfCoefficients); for (int i = 0; i < mNumberOfCoefficients; i++) { readLine(buffs, stream); mCoefficients[i] = buffs.Atof(); } + // read precision + readLine(buffs,stream); + mPrecision = buffs.Atof(); + // check end_of_data record readLine(buffs, stream); if (!buffs.BeginsWith("END") || !buffs.Contains(GetName())) { @@ -380,19 +387,24 @@ void Chebyshev3DCalc::initializeRows(int nr) { if (mNumberOfColumnsAtRow) { delete[] mNumberOfColumnsAtRow; + mNumberOfColumnsAtRow = 0; } if (mColumnAtRowBeginning) { delete[] mColumnAtRowBeginning; + mColumnAtRowBeginning = 0; } if (mTemporaryCoefficients1D) { delete[] mTemporaryCoefficients1D; + mTemporaryCoefficients1D = 0; } mNumberOfRows = nr; - mNumberOfColumnsAtRow = new UShort_t[mNumberOfRows]; - mTemporaryCoefficients1D = new Float_t[mNumberOfRows]; - mColumnAtRowBeginning = new UShort_t[mNumberOfRows]; - for (int i = mNumberOfRows; i--;) { - mNumberOfColumnsAtRow[i] = mColumnAtRowBeginning[i] = 0; + if (mNumberOfRows) { + mNumberOfColumnsAtRow = new UShort_t[mNumberOfRows]; + mTemporaryCoefficients1D = new Float_t[mNumberOfRows]; + mColumnAtRowBeginning = new UShort_t[mNumberOfRows]; + for (int i = mNumberOfRows; i--;) { + mNumberOfColumnsAtRow[i] = mColumnAtRowBeginning[i] = 0; + } } } @@ -401,23 +413,28 @@ void Chebyshev3DCalc::initializeColumns(int nc) mNumberOfColumns = nc; if (mTemporaryCoefficients2D) { delete[] mTemporaryCoefficients2D; + mTemporaryCoefficients2D = 0; } - mTemporaryCoefficients2D = new Float_t[mNumberOfColumns]; + if (mNumberOfColumns) mTemporaryCoefficients2D = new Float_t[mNumberOfColumns]; } void Chebyshev3DCalc::initializeElementBound2D(int ne) { if (mCoefficientBound2D0) { delete[] mCoefficientBound2D0; + mCoefficientBound2D0 = 0; } if (mCoefficientBound2D1) { delete[] mCoefficientBound2D1; + mCoefficientBound2D1 = 0; } mNumberOfElementsBound2D = ne; - mCoefficientBound2D0 = new UShort_t[mNumberOfElementsBound2D]; - mCoefficientBound2D1 = new UShort_t[mNumberOfElementsBound2D]; - for (int i = mNumberOfElementsBound2D; i--;) { - mCoefficientBound2D0[i] = mCoefficientBound2D1[i] = 0; + if (mNumberOfElementsBound2D) { + mCoefficientBound2D0 = new UShort_t[mNumberOfElementsBound2D]; + mCoefficientBound2D1 = new UShort_t[mNumberOfElementsBound2D]; + for (int i = mNumberOfElementsBound2D; i--;) { + mCoefficientBound2D0[i] = mCoefficientBound2D1[i] = 0; + } } } @@ -425,11 +442,14 @@ void Chebyshev3DCalc::initializeCoefficients(int nc) { if (mCoefficients) { delete[] mCoefficients; + mCoefficients = 0; } mNumberOfCoefficients = nc; - mCoefficients = new Float_t[mNumberOfCoefficients]; - for (int i = mNumberOfCoefficients; i--;) { - mCoefficients[i] = 0.0; + if (mNumberOfCoefficients) { + mCoefficients = new Float_t[mNumberOfCoefficients]; + for (int i = mNumberOfCoefficients; i--;) { + mCoefficients[i] = 0.0; + } } } diff --git a/field/Chebyshev3DCalc.h b/MathUtils/Chebyshev3DCalc.h similarity index 93% rename from field/Chebyshev3DCalc.h rename to MathUtils/Chebyshev3DCalc.h index 9c7b9419a55a1..35eea3eea7682 100644 --- a/field/Chebyshev3DCalc.h +++ b/MathUtils/Chebyshev3DCalc.h @@ -2,21 +2,21 @@ /// \brief Definition of the Cheb3DCalc class /// \author ruben.shahoyan@cern.ch 09/09/2006 -#ifndef ALICEO2_FIELD_CHEBYSHEV3DCALC_H_ -#define ALICEO2_FIELD_CHEBYSHEV3DCALC_H_ +#ifndef ALICEO2_MATHUTILS_CHEBYSHEV3DCALC_H_ +#define ALICEO2_MATHUTILS_CHEBYSHEV3DCALC_H_ #include class TSystem; // To decrease the compilable code size comment this define. This will exclude the routines // used for the calculation and saving of the coefficients. -//#define _INC_CREATION_ALICHEB3D_ +#define _INC_CREATION_Chebyshev3D_ // When _BRING_TO_BOUNDARY_ is defined, the point outside of the fitted folume is assumed to be on the surface // #define _BRING_TO_BOUNDARY_ namespace AliceO2 { -namespace Field { +namespace MathUtils { class Chebyshev3DCalc : public TNamed { public: @@ -52,7 +52,7 @@ class Chebyshev3DCalc : public TNamed { /// VERY IMPORTANT: par must contain the function arguments ALREADY MAPPED to [-1:1] interval Float_t evaluateDerivative2(int dim1, int dim2, const Float_t* par) const; -#ifdef _INC_CREATION_ALICHEB3D_ +#ifdef _INC_CREATION_Chebyshev3D_ /// Writes coefficients data to output text file, optionally appending on the end of existing file void saveData(const char* outfile, Bool_t append = kFALSE) const; @@ -95,11 +95,22 @@ class Chebyshev3DCalc : public TNamed { return mNumberOfColumnsAtRow; } - UShort_t* GetColAtRowBg() const + UShort_t* getColAtRowBg() const { return mColumnAtRowBeginning; } + Float_t getPrecision() const + { + return mPrecision; + } + + /// Sets requested precision + void setPrecision(Float_t prc=1e-6) + { + mPrecision = prc; + } + /// Sets maximum number of significant coefficients for given row/column of coefficients 3D matrix void initializeElementBound2D(int ne); @@ -143,6 +154,7 @@ class Chebyshev3DCalc : public TNamed { Int_t mNumberOfRows; ///< number of significant rows in the 3D coeffs matrix Int_t mNumberOfColumns; ///< max number of significant cols in the 3D coeffs matrix Int_t mNumberOfElementsBound2D; ///< number of elements (mNumberOfRows*mNumberOfColumns) to store for the 2D boundary + Float_t mPrecision; ///< requested precision /// of significant coeffs UShort_t* mNumberOfColumnsAtRow; //[mNumberOfRows] number of sighificant columns (2nd dim) at each row of 3D coefs matrix @@ -158,7 +170,7 @@ class Chebyshev3DCalc : public TNamed { Float_t* mTemporaryCoefficients2D; //[mNumberOfColumns] temp. coeffs for 2d summation Float_t* mTemporaryCoefficients1D; //[mNumberOfRows] temp. coeffs for 1d summation - ClassDef(AliceO2::Field::Chebyshev3DCalc, 2) // Class for interpolation of 3D->1 function by Chebyshev parametrization + ClassDef(AliceO2::MathUtils::Chebyshev3DCalc, 2) // Class for interpolation of 3D->1 function by Chebyshev parametrization }; /// Evaluates 1D Chebyshev parameterization. x is the argument mapped to [-1:1] interval diff --git a/MathUtils/MathUtilsLinkDef.h b/MathUtils/MathUtilsLinkDef.h new file mode 100644 index 0000000000000..62394e7fa058d --- /dev/null +++ b/MathUtils/MathUtilsLinkDef.h @@ -0,0 +1,10 @@ + +#ifdef __CINT__ +#pragma link off all globals; +#pragma link off all classes; +#pragma link off all functions; + +#pragma link C++ class AliceO2::MathUtils::Chebyshev3D+; +#pragma link C++ class AliceO2::MathUtils::Chebyshev3DCalc+; + +#endif diff --git a/Resources/maps/mfchebKGI_meas.root b/Resources/maps/mfchebKGI_meas.root deleted file mode 100644 index 5e9ec40ff771d..0000000000000 Binary files a/Resources/maps/mfchebKGI_meas.root and /dev/null differ diff --git a/Resources/maps/mfchebKGI_sym.root b/Resources/maps/mfchebKGI_sym.root index 142e628c4baee..01b045f527df6 100644 Binary files a/Resources/maps/mfchebKGI_sym.root and b/Resources/maps/mfchebKGI_sym.root differ diff --git a/field/CMakeLists.txt b/field/CMakeLists.txt index 7275e04ac966d..907171d46a246 100644 --- a/field/CMakeLists.txt +++ b/field/CMakeLists.txt @@ -3,6 +3,7 @@ # The extension is already found. Any number of sources could be listed here. set(INCLUDE_DIRECTORIES +${CMAKE_SOURCE_DIR} ${CMAKE_SOURCE_DIR}/field ${CMAKE_SOURCE_DIR}/header ${BASE_INCLUDE_DIRECTORIES} @@ -12,8 +13,6 @@ ${ROOT_INCLUDE_DIR} include_directories( ${INCLUDE_DIRECTORIES}) set(LINK_DIRECTORIES -${CMAKE_SOURCE_DIR}/base -${CMAKE_SOURCE_DIR}/field ${FAIRROOT_LIBRARY_DIR} ${ROOT_LIBRARY_DIR} ) @@ -23,14 +22,12 @@ link_directories( ${LINK_DIRECTORIES}) set(SRCS MagneticWrapperChebyshev.cxx MagneticField.cxx -Chebyshev3D.cxx -Chebyshev3DCalc.cxx ) Set(HEADERS) Set(LINKDEF fieldLinkDef.h) Set(LIBRARY_NAME Field) -Set(DEPENDENCIES Base EG Physics Cint Core) +Set(DEPENDENCIES MathUtils Cint Core) GENERATE_LIBRARY() diff --git a/field/MagneticWrapperChebyshev.cxx b/field/MagneticWrapperChebyshev.cxx index ad338015c21ff..e65922fd79718 100644 --- a/field/MagneticWrapperChebyshev.cxx +++ b/field/MagneticWrapperChebyshev.cxx @@ -9,6 +9,7 @@ #include "FairLogger.h" using namespace AliceO2::Field; +using namespace AliceO2::MathUtils; ClassImp(MagneticWrapperChebyshev) @@ -788,7 +789,7 @@ void MagneticWrapperChebyshev::getTPCRatIntegralCylindrical(const Double_t* rphi return; } -#ifdef _INC_CREATION_ALICHEB3D_ +#ifdef _INC_CREATION_Chebyshev3D_ void MagneticWrapperChebyshev::loadData(const char* inpfile) { @@ -847,7 +848,7 @@ void MagneticWrapperChebyshev::loadData(const char* inpfile) for (int ip = 0; ip < nparTPCInt; ip++) { Chebyshev3D* cheb = new Chebyshev3D(); cheb->loadData(stream); - AddParamTPCInt(cheb); + addParameterTPCIntegral(cheb); } Chebyshev3DCalc::readLine(buffs, stream); @@ -871,7 +872,7 @@ void MagneticWrapperChebyshev::loadData(const char* inpfile) for (int ip = 0; ip < nparTPCRatInt; ip++) { Chebyshev3D* cheb = new Chebyshev3D(); cheb->loadData(stream); - AddParamTPCRatInt(cheb); + addParameterTPCRatIntegral(cheb); } Chebyshev3DCalc::readLine(buffs, stream); @@ -907,8 +908,8 @@ void MagneticWrapperChebyshev::loadData(const char* inpfile) Chebyshev3DCalc::readLine(buffs, stream); - if (!buffs.BeginsWith("END DIPOLE") || !buffs.Contains(GetName())) { - Error("LoadData", "Expected: \"END DIPOLE\", found \"%s\"\nStop\n", buffs.Data()); + if (!buffs.BeginsWith("END ") && !buffs.Contains(GetName())) { + Error("LoadData", "Expected: \"END %s\", found \"%s\"\nStop\n", GetName(), buffs.Data()); exit(1); } @@ -958,7 +959,7 @@ void MagneticWrapperChebyshev::buildTableTPCRatIntegral() #endif -#ifdef _INC_CREATION_ALICHEB3D_ +#ifdef _INC_CREATION_Chebyshev3D_ MagneticWrapperChebyshev::MagneticWrapperChebyshev(const char* inputFile) : mNumberOfParameterizationSolenoid(0), mNumberOfDistinctZSegmentsSolenoid(0), @@ -1063,7 +1064,7 @@ void MagneticWrapperChebyshev::addParameterTPCRatIntegral(const Chebyshev3D* par } } -void MagneticWrapperChebyshev::AddParamDipole(const Chebyshev3D* param) +void MagneticWrapperChebyshev::addParameterDipole(const Chebyshev3D* param) { if (!mParameterizationDipole) { mParameterizationDipole = new TObjArray(); diff --git a/field/MagneticWrapperChebyshev.h b/field/MagneticWrapperChebyshev.h index d7aa2e8c12c5a..835ac870306fe 100644 --- a/field/MagneticWrapperChebyshev.h +++ b/field/MagneticWrapperChebyshev.h @@ -8,7 +8,7 @@ #include #include #include -#include "Chebyshev3D.h" +#include "MathUtils/Chebyshev3D.h" class TSystem; class TArrayF; @@ -171,24 +171,24 @@ class MagneticWrapperChebyshev : public TNamed { return mMaxRadiusTPCRat; } - Chebyshev3D* getParameterSolenoid(Int_t ipar) const + AliceO2::MathUtils::Chebyshev3D* getParameterSolenoid(Int_t ipar) const { - return (Chebyshev3D*)mParameterizationSolenoid->UncheckedAt(ipar); + return (AliceO2::MathUtils::Chebyshev3D*)mParameterizationSolenoid->UncheckedAt(ipar); } - Chebyshev3D* getParameterTPCRatIntegral(Int_t ipar) const + AliceO2::MathUtils::Chebyshev3D* getParameterTPCRatIntegral(Int_t ipar) const { - return (Chebyshev3D*)mParameterizationTPCRat->UncheckedAt(ipar); + return (AliceO2::MathUtils::Chebyshev3D*)mParameterizationTPCRat->UncheckedAt(ipar); } - Chebyshev3D* getParameterTPCIntegral(Int_t ipar) const + AliceO2::MathUtils::Chebyshev3D* getParameterTPCIntegral(Int_t ipar) const { - return (Chebyshev3D*)mParameterizationTPC->UncheckedAt(ipar); + return (AliceO2::MathUtils::Chebyshev3D*)mParameterizationTPC->UncheckedAt(ipar); } - Chebyshev3D* getParameterDipole(Int_t ipar) const + AliceO2::MathUtils::Chebyshev3D* getParameterDipole(Int_t ipar) const { - return (Chebyshev3D*)mParameterizationDipole->UncheckedAt(ipar); + return (AliceO2::MathUtils::Chebyshev3D*)mParameterizationDipole->UncheckedAt(ipar); } /// Prints info @@ -237,7 +237,7 @@ class MagneticWrapperChebyshev : public TNamed { static void cartesianToCylindrical(const Double_t* xyz, Double_t* rphiz); static void cylindricalToCartesian(const Double_t* rphiz, Double_t* xyz); -#ifdef _INC_CREATION_ALICHEB3D_ // see Cheb3D.h for explanation +#ifdef _INC_CREATION_Chebyshev3D_ // see Cheb3D.h for explanation /// Reads coefficients data from the text file void loadData(const char* inpfile); @@ -254,18 +254,18 @@ class MagneticWrapperChebyshev : public TNamed { /// Adds new parameterization piece for Solenoid /// NOTE: pieces must be added strictly in increasing R then increasing Z order - void addParameterSolenoid(const Chebyshev3D* param); + void addParameterSolenoid(const AliceO2::MathUtils::Chebyshev3D* param); // Adds new parameterization piece for TPCIntegral // NOTE: pieces must be added strictly in increasing R then increasing Z order - void addParameterTPCIntegral(const Chebyshev3D* param); + void addParameterTPCIntegral(const AliceO2::MathUtils::Chebyshev3D* param); /// Adds new parameterization piece for TPCRatInt // NOTE: pieces must be added strictly in increasing R then increasing Z order - void addParameterTPCRatIntegral(const Chebyshev3D* param); + void addParameterTPCRatIntegral(const AliceO2::MathUtils::Chebyshev3D* param); /// Adds new parameterization piece for Dipole - void addParameterDipole(const Chebyshev3D* param); + void addParameterDipole(const AliceO2::MathUtils::Chebyshev3D* param); /// Builds lookup table for dipole void buildTable(Int_t npar, TObjArray* parArr, Int_t& nZSeg, Int_t& nYSeg, Int_t& nXSeg, Float_t& minZ, Float_t& maxZ, @@ -395,7 +395,7 @@ class MagneticWrapperChebyshev : public TNamed { Float_t mMaxDipoleZ; ///< Max Z of Dipole parameterization TObjArray* mParameterizationDipole; ///< Parameterization pieces for Dipole field - FairLogger* mLogger; + FairLogger* mLogger; //! ClassDef(AliceO2::Field::MagneticWrapperChebyshev, 2) // Wrapper class for the set of Chebishev parameterizations of Alice mag.field }; diff --git a/field/fieldLinkDef.h b/field/fieldLinkDef.h index c340fb3d5f496..a3bc529576e02 100644 --- a/field/fieldLinkDef.h +++ b/field/fieldLinkDef.h @@ -4,89 +4,7 @@ #pragma link off all classes; #pragma link off all functions; -#pragma link C++ class AliceO2::Field::Chebyshev3D+; #pragma link C++ class AliceO2::Field::MagneticField+; #pragma link C++ class AliceO2::Field::MagneticWrapperChebyshev+; -#pragma link C++ class AliceO2::Field::Chebyshev3DCalc+; - -//map the old to the new variables: -//AliCheb3D-------OLD -//Int_t fDimOut; // dimension of the ouput array -//Float_t fPrec; // requested precision -//Float_t fBMin[3]; // min boundaries in each dimension -//Float_t fBMax[3]; // max boundaries in each dimension -//Float_t fBScale[3]; // scale for boundary mapping to [-1:1] interval -//Float_t fBOffset[3]; // offset for boundary mapping to [-1:1] interval -//TObjArray fChebCalc; // Chebyshev parameterization for each output dimension - -//Chebyshev3D-------New -//Int_t mOutputArrayDimension; ///< dimension of the ouput array -//Float_t mPrecision; ///< requested precision -//Float_t mMinBoundaries[3]; ///< min boundaries in each dimension -//Float_t mMaxBoundaries[3]; ///< max boundaries in each dimension -//Float_t mBoundaryMappingScale[3]; ///< scale for boundary mapping to [-1:1] interval -//Float_t mBoundaryMappingOffset[3]; ///< offset for boundary mapping to [-1:1] interval -//TObjArray mChebyshevParameter; ///< Chebyshev parameterization for each output dimension - -#pragma read sourceClass="AliCheb3D" \ - targetClass="AliceO2::Field::Chebyshev3D" \ - version="[1-]" \ - source="Int_t fDimOut" \ - target="mOutputArrayDimension" \ - code="{ mOutputArrayDimension = onfile.fDimOut; }" - -#pragma read sourceClass="AliCheb3D" \ - targetClass="AliceO2::Field::Chebyshev3D" \ - version="[1-]" \ - source="Float_t fPrec" \ - target="mPrecision" \ - code="{ mPrecision = onfile.fPrec; }" - - -#pragma read sourceClass="AliCheb3D" \ - targetClass="AliceO2::Field::Chebyshev3D" \ - version="[1-]" \ - source="Float_t fBMin[3]" \ - target="mMinBoundaries" \ - code="{ mMinBoundaries[0] = onfile.fBMin[0];mMinBoundaries[1] = onfile.fBMin[1];mMinBoundaries[2] = onfile.fBMin[2]; }" - -#pragma read sourceClass="AliCheb3D" \ - targetClass="AliceO2::Field::Chebyshev3D" \ - version="[1-]" \ - source="Float_t fBMax[3]" \ - target="mMaxBoundaries" \ - code="{ mMaxBoundaries[0] = onfile.fBMax[0];mMaxBoundaries[1] = onfile.fBMax[1];mMaxBoundaries[2] = onfile.fBMax[2]; }" - -#pragma read sourceClass="AliCheb3D" \ - targetClass="AliceO2::Field::Chebyshev3D" \ - version="[1-]" \ - source="Float_t fBScale[3]" \ - target="mBoundaryMappingScale" \ - code="{ mBoundaryMappingScale[0] = onfile.fBScale[0];mBoundaryMappingScale[1] = onfile.fBScale[1];mBoundaryMappingScale[2] = onfile.fBScale[2]; }" - -#pragma read sourceClass="AliCheb3D" \ - targetClass="AliceO2::Field::Chebyshev3D" \ - version="[1-]" \ - source="Float_t fBOffset[3]" \ - target="mBoundaryMappingScale" \ - code="{ mBoundaryMappingScale[0] = onfile.fBOffset[0];mBoundaryMappingScale[1] = onfile.fBOffset[1];mBoundaryMappingScale[2] = onfile.fBOffset[2]; }" - -#pragma read sourceClass="AliCheb3D" \ - targetClass="AliceO2::Field::Chebyshev3D" \ - version="[1-]" \ - source="TObjArray fChebCalc" \ - target="mChebyshevParameter" \ - code="{ mChebyshevParameter = onfile.fChebCalc;}" - - - - - - - - - - - #endif