Jpp test-rotations-old
the software that should make you happy
Loading...
Searching...
No Matches
JACOUSTICS::JKatoomba< JGandalf > Struct Reference

Template specialisation of fit function of acoustic model based on JGandalf minimiser. More...

#include <JKatoomba_t.hh>

Inheritance diagram for JACOUSTICS::JKatoomba< JGandalf >:
JACOUSTICS::JKatoomba< JMinimiser_t >

Public Types

typedef double result_type
 
typedef std::map< JLocation, std::map< JEmitter, std::vector< JHit > > > data_type
 Type definition internal data structure.
 

Public Member Functions

 JKatoomba (const JGeometry &geometry, const JSoundVelocity &velocity, const int option)
 Constructor The option corresponds to the use of fit parameters in the model of the detector geometry.
 
template<class T >
result_type operator() (T __begin, T __end)
 Fit.
 

Public Attributes

double lambda
 
JModel value
 
int numberOfIterations
 
JMATH::JMatrixNS V
 

Static Public Attributes

static int debug
 debug level
 
static int MAXIMUM_ITERATIONS
 maximal number of iterations
 
static double EPSILON
 maximal distance to minimum
 
static double LAMBDA_MIN
 minimal value control parameter
 
static double LAMBDA_MAX
 maximal value control parameter
 
static double LAMBDA_UP
 multiplication factor control parameter
 
static double LAMBDA_DOWN
 multiplication factor control parameter
 
static double PIVOT
 minimal value diagonal element of matrix
 

Private Member Functions

void evaluate (const data_type &data)
 Evaluation of fit.
 

Private Attributes

JMATH::JVectorND Y
 
result_type successor
 
JModel previous
 
std::vector< double > h
 

Detailed Description

Template specialisation of fit function of acoustic model based on JGandalf minimiser.

Definition at line 677 of file JKatoomba_t.hh.

Member Typedef Documentation

◆ result_type

double JACOUSTICS::JKatoomba< JGandalf >::result_type

Definition at line 680 of file JKatoomba_t.hh.

◆ data_type

Type definition internal data structure.

Definition at line 685 of file JKatoomba_t.hh.

Constructor & Destructor Documentation

◆ JKatoomba()

JACOUSTICS::JKatoomba< JGandalf >::JKatoomba ( const JGeometry & geometry,
const JSoundVelocity & velocity,
const int option )
inline

Constructor The option corresponds to the use of fit parameters in the model of the detector geometry.


A negative implies that all strings in the detector use common fit parameters.

Parameters
geometrydetector geometry
velocitysound velocity
optionoption

Definition at line 697 of file JKatoomba_t.hh.

699 :
700 JKatoomba<>(geometry, velocity, option)
701 {};
Template definition of fit function of acoustic model.

Member Function Documentation

◆ operator()()

template<class T >
result_type JACOUSTICS::JKatoomba< JGandalf >::operator() ( T __begin,
T __end )
inline

Fit.

Parameters
__beginbegin of hits
__endend of hits
Returns
chi2 and gradient

Definition at line 712 of file JKatoomba_t.hh.

713 {
714 using namespace std;
715 using namespace JPP;
716
717 value.setOption(this->option);
718
719 const int N = value.getN();
720
721 V.resize(N);
722 Y.resize(N);
723 h.resize(N);
724
726
727 for (T hit = __begin; hit != __end; ++hit) {
728 data[hit->getLocation()][hit->getEmitter()].push_back(*hit);
729 }
730
732
733 result_type precessor = numeric_limits<double>::max();
734
736
737 DEBUG("step: " << numberOfIterations << endl);
738
739 evaluate(data);
740
741 DEBUG("lambda: " << FIXED(12,5) << lambda << endl);
742 DEBUG("chi2: " << FIXED(12,5) << successor << endl);
743
744 if (successor < precessor) {
745
746 if (numberOfIterations != 0) {
747
748 if (fabs(precessor - successor) < EPSILON*fabs(precessor)) {
749 return successor;
750 }
751
752 if (lambda > LAMBDA_MIN) {
754 }
755 }
756
757 precessor = successor;
758 previous = value;
759
760 } else {
761
762 value = previous;
763 lambda *= LAMBDA_UP;
764
765 if (lambda > LAMBDA_MAX) {
766 return precessor; // no improvement found
767 }
768
769 evaluate(data);
770 }
771
772 // force definite positiveness
773
774 for (int i = 0; i != N; ++i) {
775
776 if (V(i,i) < PIVOT) {
777 V(i,i) = PIVOT;
778 }
779
780 h[i] = 1.0 / sqrt(V(i,i));
781 }
782
783 // normalisation
784
785 for (int i = 0; i != N; ++i) {
786 for (int j = 0; j != i; ++j) {
787 V(j,i) *= h[i] * h[j];
788 V(i,j) = V(j,i);
789 }
790 }
791
792 for (int i = 0; i != N; ++i) {
793 V(i,i) = 1.0 + lambda;
794 }
795
796
797 // solve A x = b
798
799 for (int col = 0; col != N; ++col) {
800 Y[col] *= h[col];
801 }
802
803 try {
804 V.solve(Y);
805 }
806 catch (const exception& error) {
807
808 ERROR("JGandalf: " << error.what() << endl << V << endl);
809
810 break;
811 }
812
813 // update value
814
815 for (int row = 0; row != N; ++row) {
816
817 DEBUG("u[" << noshowpos << setw(3) << row << "] = " << showpos << FIXED(20,5) << value[row]);
818
819 value[row] -= h[row] * Y[row];
820
821 DEBUG(" -> " << FIXED(20,10) << value[row] << noshowpos << endl);
822 }
823 }
824
825 return precessor;
826 }
#define DEBUG(A)
Message macros.
Definition JMessage.hh:62
#define ERROR(A)
Definition JMessage.hh:66
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
int j
Definition JPolint.hh:801
Auxiliary data structure for floating point format specification.
Definition JManip.hh:448
static double LAMBDA_DOWN
multiplication factor control parameter
static double PIVOT
minimal value diagonal element of matrix
static double LAMBDA_UP
multiplication factor control parameter
static int MAXIMUM_ITERATIONS
maximal number of iterations
static double LAMBDA_MAX
maximal value control parameter
static double EPSILON
maximal distance to minimum
void evaluate(const data_type &data)
Evaluation of fit.
static double LAMBDA_MIN
minimal value control parameter
void setOption(const int option)
Set fit option.
size_t getN() const
Get number of fit parameters.
Data structure for measured coincidence rates of all pairs of PMTs in optical module.
Definition JFitK40.hh:103
void resize(const size_t size)
Resize matrix.
Definition JMatrixND.hh:446
void solve(JVectorND_t &u)
Get solution of equation A x = b.
Definition JMatrixNS.hh:308

◆ evaluate()

void JACOUSTICS::JKatoomba< JGandalf >::evaluate ( const data_type & data)
inlineprivate

Evaluation of fit.

Parameters
datadata

Definition at line 848 of file JKatoomba_t.hh.

849 {
850 using namespace std;
851 using namespace JPP;
852
853 successor = 0.0;
854
855 V.reset();
856 Y.reset();
857
858 for (data_type::const_iterator p = data.begin(); p != data.end(); ++p) {
859
860 const JGEOMETRY::JString& string = geometry [p->first.getString()];
861 const JMODEL ::JString& parameters = value.string[p->first.getString()];
862 const JPosition3D position = string.getPosition(parameters, p->first.getFloor());
863
864 for (data_type::mapped_type::const_iterator emitter = p->second.begin(); emitter != p->second.end(); ++emitter) {
865
866 const double D = emitter->first.getDistance(position);
867 const double Vi = velocity.getInverseVelocity(D, emitter->first.getZ(), position.getZ());
868
869 const H_t H0(1.0, string.getGradient(parameters, emitter->first.getPosition(), p->first.getFloor()) * Vi);
870
871 for (data_type::mapped_type::mapped_type::const_iterator hit = emitter->second.begin(); hit != emitter->second.end(); ++hit) {
872
873 const double toa_s = value.emission[hit->getEKey()].t1 + D * Vi;
874
875 const double u = (toa_s - hit->getValue()) / hit->getSigma();
876 const double W = sqrt(hit->getWeight());
877
878 successor += (W*W) * estimator->getRho(u);
879
880 const H_t H = H0 * (W * estimator->getPsi(u) / hit->getSigma());
881
882 I_t i;
883
884 i.t1 = value.getIndex(hit->getEKey(), &H_t::t1);
885 i.tx = value.getIndex(hit->getString(), &H_t::tx);
886 i.ty = value.getIndex(hit->getString(), &H_t::ty);
887 i.tx2 = value.getIndex(hit->getString(), &H_t::tx2);
888 i.ty2 = value.getIndex(hit->getString(), &H_t::ty2);
889 i.vs = value.getIndex(hit->getString(), &H_t::vs);
890
891 V(i.t1, i.t1) += H.t1 * H.t1;
892
893 Y[i.t1] += W * H.t1;
894
895 if (hit->getFloor() != 0) {
896
897 switch (this->option) {
898
900 V(i.t1, i.vs) += H.t1 * H.vs; V(i.tx, i.vs) += H.tx * H.vs; V(i.ty, i.vs) += H.ty * H.vs; V(i.tx2, i.vs) += H.tx2 * H.vs; V(i.ty2, i.vs) += H.ty2 * H.vs;
901
902 V(i.vs, i.t1) = V(i.t1, i.vs);
903 V(i.vs, i.tx) = V(i.tx, i.vs);
904 V(i.vs, i.ty) = V(i.ty, i.vs);
905 V(i.vs, i.tx2) = V(i.tx2, i.vs);
906 V(i.vs, i.ty2) = V(i.ty2, i.vs);
907
908 V(i.vs, i.vs) += H.vs * H.vs;
909
910 Y[i.vs] += W * H.vs;
911
913 V(i.t1, i.tx2) += H.t1 * H.tx2; V(i.tx, i.tx2) += H.tx * H.tx2; V(i.ty, i.tx2) += H.ty * H.tx2;
914
915 V(i.tx2, i.t1) = V(i.t1, i.tx2);
916 V(i.tx2, i.tx) = V(i.tx, i.tx2);
917 V(i.tx2, i.ty) = V(i.ty, i.tx2);
918
919 V(i.t1, i.ty2) += H.t1 * H.ty2; V(i.tx, i.ty2) += H.tx * H.ty2; V(i.ty, i.ty2) += H.ty * H.ty2;
920
921 V(i.ty2, i.t1) = V(i.t1, i.ty2);
922 V(i.ty2, i.tx) = V(i.tx, i.ty2);
923 V(i.ty2, i.ty) = V(i.ty, i.ty2);
924
925 V(i.tx2, i.tx2) += H.tx2 * H.tx2; V(i.tx2, i.ty2) += H.tx2 * H.ty2;
926 V(i.ty2, i.tx2) = V(i.tx2, i.ty2); V(i.ty2, i.ty2) += H.ty2 * H.ty2;
927
928 Y[i.tx2] += W * H.tx2;
929 Y[i.ty2] += W * H.ty2;
930
932 V(i.t1, i.tx) += H.t1 * H.tx; V(i.t1, i.ty) += H.t1 * H.ty;
933 V(i.tx, i.t1) = V(i.t1, i.tx); V(i.ty, i.t1) = V(i.t1, i.ty);
934
935 V(i.tx, i.tx) += H.tx * H.tx; V(i.tx, i.ty) += H.tx * H.ty;
936 V(i.ty, i.tx) = V(i.tx, i.ty); V(i.ty, i.ty) += H.ty * H.ty;
937
938 Y[i.tx] += W * H.tx;
939 Y[i.ty] += W * H.ty;
940 break;
941
942 default:
943 break;
944 }
945 }
946 }
947 }
948 }
949 }
Data structure for position in three dimensions.
const JPosition3D & getPosition() const
Get position.
double getZ() const
Get z position.
Definition JVector3D.hh:115
@ FIT_EMITTERS_AND_STRINGS_1st_ORDER_t
fit times of emission of emitters and tilt angles of strings
@ FIT_EMITTERS_AND_STRINGS_2nd_ORDER_t
fit times of emission of emitters and tilt angles of strings with second order correction
@ FIT_EMITTERS_AND_STRINGS_2nd_ORDER_AND_STRETCHING_t
fit times of emission of emitters and tilt angles of strings with second order correction and stretch...
static const double H
Planck constant [eV s].
JACOUSTICS::JModel::emission_type emission
size_t getIndex(int id, double JString::*p) const
Get index of fit parameter for given string.
JACOUSTICS::JModel::string_type string
JMatrixND & reset()
Set matrix to the null matrix.
Definition JMatrixND.hh:459
void reset()
Reset.
Definition JVectorND.hh:45

Member Data Documentation

◆ debug

int JACOUSTICS::JKatoomba< JGandalf >::debug
static

debug level

Definition at line 828 of file JKatoomba_t.hh.

◆ MAXIMUM_ITERATIONS

int JACOUSTICS::JKatoomba< JGandalf >::MAXIMUM_ITERATIONS
static

maximal number of iterations

Definition at line 829 of file JKatoomba_t.hh.

◆ EPSILON

double JACOUSTICS::JKatoomba< JGandalf >::EPSILON
static

maximal distance to minimum

Definition at line 830 of file JKatoomba_t.hh.

◆ LAMBDA_MIN

double JACOUSTICS::JKatoomba< JGandalf >::LAMBDA_MIN
static

minimal value control parameter

Definition at line 831 of file JKatoomba_t.hh.

◆ LAMBDA_MAX

double JACOUSTICS::JKatoomba< JGandalf >::LAMBDA_MAX
static

maximal value control parameter

Definition at line 832 of file JKatoomba_t.hh.

◆ LAMBDA_UP

double JACOUSTICS::JKatoomba< JGandalf >::LAMBDA_UP
static

multiplication factor control parameter

Definition at line 833 of file JKatoomba_t.hh.

◆ LAMBDA_DOWN

double JACOUSTICS::JKatoomba< JGandalf >::LAMBDA_DOWN
static

multiplication factor control parameter

Definition at line 834 of file JKatoomba_t.hh.

◆ PIVOT

double JACOUSTICS::JKatoomba< JGandalf >::PIVOT
static

minimal value diagonal element of matrix

Definition at line 835 of file JKatoomba_t.hh.

◆ lambda

double JACOUSTICS::JKatoomba< JGandalf >::lambda

Definition at line 837 of file JKatoomba_t.hh.

◆ value

Definition at line 838 of file JKatoomba_t.hh.

◆ numberOfIterations

int JACOUSTICS::JKatoomba< JGandalf >::numberOfIterations

Definition at line 839 of file JKatoomba_t.hh.

◆ V

◆ Y

Definition at line 952 of file JKatoomba_t.hh.

◆ successor

Definition at line 953 of file JKatoomba_t.hh.

◆ previous

JModel JACOUSTICS::JKatoomba< JGandalf >::previous
private

Definition at line 954 of file JKatoomba_t.hh.

◆ h

Definition at line 955 of file JKatoomba_t.hh.


The documentation for this struct was generated from the following file: