ATLAS Offline Software
Loading...
Searching...
No Matches
Calibrator.cxx File Reference
#include "Calibrator.h"
#include "CxxUtils/checker_macros.h"
#include <TH1F.h>
#include <TH2F.h>
#include <TH1D.h>
#include <TDirectory.h>
#include <TF1.h>
#include <TStyle.h>
#include <TNtuple.h>
#include <TGraphErrors.h>
#include <iostream>
#include <iomanip>
#include <fstream>
#include <cmath>

Go to the source code of this file.

Functions

double pol3deg (double *x, double *par)
double pol3deg2 (double *x, double *par)
double trrel_dines (double *x, double *par)
double rtrel_dines (double *x, double *par)
float Calibrator::FitRt ATLAS_NOT_THREAD_SAFE (const std::string &key, const std::string &opt, TH2F *rtHist, TDirectory *dir)
float Calibrator::FitTimeResidual ATLAS_NOT_THREAD_SAFE (const std::string &key, TH1F *tresHist)
TDirectory *Calibrator::Calibrate ATLAS_NOT_THREAD_SAFE (TDirectory *dir, std::string key, const std::string &opt, caldata *caldata_above)

Function Documentation

◆ ATLAS_NOT_THREAD_SAFE() [1/3]

float Calibrator::FitRt ATLAS_NOT_THREAD_SAFE ( const std::string & key,
const std::string & opt,
TH2F * rtHist,
TDirectory * dir )

Definition at line 555 of file Calibrator.cxx.

555 { // Global gStyle is used.
556
557 float rtpars[4];
558
559 //create r-m_t and m_t-r graphs
560 RtGraph* rtg = m_rtbinning.find('t')==std::string::npos ? new RtGraph(rtHist,1,std::string("abs(rtrack)").data(),!bequiet,dir) : new RtGraph(rtHist,0,std::string("t-t0").data(),!bequiet,dir);
561
562 TF1 dtfitfunc("dtfitfunc","gaus(0)");
563
564 gStyle->SetOptFit(1111);
565 gStyle->SetPalette(1);
566
567 // select type of R-m_t relation
568 TF1 *trfunc,*rtfunc,*oldrtfunc;
569 if (m_isdines){
570 trfunc = new TF1("trfunc",trrel_dines,0,3,4);
571 rtfunc = new TF1("rtfunc",rtrel_dines,rtg->mintime,200,4);
572 oldrtfunc = new TF1("oldrtfunc",rtrel_dines,rtg->mintime,200,4);
573 trfunc->SetParameters(0.3, 1.0, 0.05, -3.);
574 double rmin0Limits[2]={0.0,2.0};
575 double rhoLimits[2]={0.0,10.};
576 double vLimits[2]={0.0,1.0};
577 trfunc->SetParLimits(0,rmin0Limits[0],rmin0Limits[1]);
578 trfunc->SetParLimits(1,rhoLimits[0],rhoLimits[1]);
579 trfunc->SetParLimits(2,vLimits[0],vLimits[1]);
580 rtg->trgr->SetTitle("Dines' R-t sqrt(..)");
581 }
582 else {
583 trfunc = new TF1("trfunc",pol3deg,0,3,4);
584 rtfunc = new TF1("rtfunc",pol3deg,rtg->mintime,200,4);
585 oldrtfunc = new TF1("oldrtfunc",pol3deg,rtg->mintime,200,4);
586 trfunc->SetParameters(0.0, 0.0, 0.0, 0.0);
587 rtg->trgr->SetTitle("Polynomial R-t");
588 }
589
590 // m_t-r relation
591 trfunc->SetRange(0,2);
592 if (opt.find('3')==std::string::npos) trfunc->FixParameter(3,0);
593 trfunc->SetLineColor(4) ;
594 rtg->trgr->Fit(trfunc,"QR"); // always fit the m_t-r relation
595 if (!bequiet) rtg->trgr->Write();
596
597
598 // r-m_t relation
599 rtfunc->SetRange(0,45);
600 if (opt.find('3')==std::string::npos) rtfunc->FixParameter(3,0);
601 rtfunc->SetLineColor(4) ;
602 if (m_isdines) {
603 rtfunc->SetParameters(trfunc->GetParameter(0),trfunc->GetParameter(1),trfunc->GetParameter(2),trfunc->GetParameter(3));
604
605 }
606 else { // else do a fit
607 rtfunc->SetParameters(0.000000e+00, 6.269950e-02, -3.370054e-04, -1.244642e-07);
608 if(rtg->rtgr != nullptr) {
609 rtg->rtgr->Fit(rtfunc,"QR"); //fit Rt first time
610 for (int ipnt=0; ipnt<rtg->rtgr->GetN(); ipnt++){ //calculate m_t-error
611 float deriv = rtfunc->Derivative(rtg->rtgr->GetX()[ipnt]);
612 if(fabs(deriv)<1.0e-24) deriv=0.05;
613 rtg->rtgr->SetPointError(ipnt , rtg->rtgr->GetEY()[ipnt]/deriv , rtg->rtgr->GetEY()[ipnt]);
614 }
615 rtg->rtgr->Fit(rtfunc,"R"); //fit again
616 }
617
618 }
619 if (!bequiet and (rtg->rtgr != nullptr)) rtg->rtgr->Write();
620
621 // old r-m_t relation
622 oldrtfunc->SetRange(0,45);
623 oldrtfunc->SetLineColor(1) ;
624 oldrtfunc->SetLineStyle(2) ;
625 oldrtfunc->SetLineWidth(1) ;
626 oldrtfunc->SetParameters(data[key].oldrtpar[0],data[key].oldrtpar[1],data[key].oldrtpar[2],data[key].oldrtpar[3]);
627
628 if (!bequiet) oldrtfunc->Write();
629
630 rtpars[0] = rtfunc->GetParameter(0);
631 rtpars[1] = rtfunc->GetParameter(1);
632 rtpars[2] = rtfunc->GetParameter(2);
633 rtpars[3] = rtfunc->GetParameter(3);
634
635 //get the t0 from the tr relation
636 float tdrift = 20;
637
638 if (!m_isdines) {
639 float rdrift = 0.0;
640 int ntries = 0;
641 float precision = 0.0001;
642 int maxtries = 500;
643 float drdt;
644 float driftradius = rtpars[0]+tdrift*(rtpars[1]+tdrift*(rtpars[2]));
645 float residual = std::abs(rdrift) - driftradius;
646 while (std::abs(residual) > precision) {
647
648 drdt = rtpars[1]+tdrift*(2*rtpars[2]);
649 if(fabs(drdt)<1.0e-24) drdt=0.05;
650 tdrift = tdrift + residual/drdt;
651
652 driftradius = rtpars[0]+tdrift*(rtpars[1]+tdrift*(rtpars[2]));
653 residual = rdrift - driftradius;
654
655 ntries = ntries + 1;
656 if (ntries>maxtries){
657 break;
658 }
659 }
660 }
661
662 if (opt.find('0')==std::string::npos) {
663 if (m_isdines){
664 data[key].rtpar[0] = rtpars[0];
665 data[key].rtpar[1] = rtpars[1];
666 data[key].rtpar[2] = rtpars[2];
667 data[key].rtpar[3] = 0;
668 } else {
669 data[key].rtpar[0] = 0;
670 data[key].rtpar[1] = rtpars[1] + 2*tdrift*rtpars[2] + 3*tdrift*tdrift*rtpars[3];
671 data[key].rtpar[2] = rtpars[2] + 3*tdrift*rtpars[3];
672 data[key].rtpar[3] = rtpars[3];
673 }
674 }
675 else {
676 data[key].rtpar[0] = rtpars[0];
677 data[key].rtpar[1] = rtpars[1];
678 data[key].rtpar[2] = rtpars[2];
679 data[key].rtpar[3] = rtpars[3];
680 }
681
682 data[key].rtgraph=rtg;
683
684 delete rtfunc;
685 delete oldrtfunc;
686 delete trfunc;
687
688 return 0.0;
689}
double trrel_dines(double *x, double *par)
double pol3deg(double *x, double *par)
double rtrel_dines(double *x, double *par)
char data[hepevt_bytes_allocation_ATLAS]
Definition HepEvt.cxx:11
A class for generating a r-t and t-r graphs by binning the 2D histograms in Calibrator::rtHists in r ...
Definition Calibrator.h:34
TGraphErrors * trgr
the t(r) graph
Definition Calibrator.h:44
TGraphErrors * rtgr
the r(t) graph
Definition Calibrator.h:45
float mintime
the minimum t-value
Definition Calibrator.h:47
STL class.
STL namespace.

◆ ATLAS_NOT_THREAD_SAFE() [2/3]

float Calibrator::FitResidual ATLAS_NOT_THREAD_SAFE ( const std::string & key,
TH1F * tresHist )

Definition at line 692 of file Calibrator.cxx.

692 { // Global gStyle is used.
693
694 float mean = tresHist->GetMean();
695 float rms = tresHist->GetRMS();
696 float val,maxval=0;
697 int maxbin=1;
698 for (int j=1; j<=tresHist->GetNbinsX()+1;j++){
699 val=tresHist->GetBinContent(j);
700 if (val>maxval){
701 maxval=val;
702 maxbin=j;
703 }
704 }
705
706
707 data[key].tresMean = 0;
708 float fitmean=tresHist->GetBinCenter(maxbin);
709 if (std::abs(fitmean)>5.0){
710 data[key].tresMean = mean;
711 data[key].t0err = rms;
712 data[key].t0fittype = 1;
713 return mean;
714 }
715
716 float fitrms=rms;
717 if (fitrms>5.0) fitrms=5.0;
718
719 tresHist->Fit("gaus","QR","",mean-rms,mean+rms);
720 if (tresHist->GetFunction("gaus")->GetParameter(1) - tresHist->GetFunction("gaus")->GetParameter(2) < m_mintres ||
721 tresHist->GetFunction("gaus")->GetParameter(1) + tresHist->GetFunction("gaus")->GetParameter(2) > m_maxtres) {
722 data[key].tresMean = mean;
723 data[key].t0err = rms;
724 data[key].t0fittype = 1;
725 return mean;
726 }
727
728 tresHist->Fit("gaus","QR","",tresHist->GetFunction("gaus")->GetParameter(1) - tresHist->GetFunction("gaus")->GetParameter(2),tresHist->GetFunction("gaus")->GetParameter(1) + tresHist->GetFunction("gaus")->GetParameter(2));
729 if (tresHist->GetFunction("gaus")->GetParameter(1) - tresHist->GetFunction("gaus")->GetParameter(2) < m_mintres ||
730 tresHist->GetFunction("gaus")->GetParameter(1) + tresHist->GetFunction("gaus")->GetParameter(2) > m_maxtres) {
731 data[key].tres = tresHist->GetFunction("gaus")->GetParameter(2);
732 data[key].tresMean = tresHist->GetFunction("gaus")->GetParameter(1);
733 data[key].t0err = tresHist->GetFunction("gaus")->GetParError(1);
734 data[key].t0fittype = 2;
735 return tresHist->GetFunction("gaus")->GetParameter(1);
736 }
737
738 tresHist->Fit("gaus","QR","",tresHist->GetFunction("gaus")->GetParameter(1) - tresHist->GetFunction("gaus")->GetParameter(2),tresHist->GetFunction("gaus")->GetParameter(1) + tresHist->GetFunction("gaus")->GetParameter(2));
739 if (tresHist->GetFunction("gaus")->GetParameter(1) - tresHist->GetFunction("gaus")->GetParameter(2) < m_mintres ||
740 tresHist->GetFunction("gaus")->GetParameter(1) + tresHist->GetFunction("gaus")->GetParameter(2) > m_maxtres) {
741 data[key].tres = tresHist->GetFunction("gaus")->GetParameter(2);
742 data[key].tresMean = tresHist->GetFunction("gaus")->GetParameter(1);
743 data[key].t0err = tresHist->GetFunction("gaus")->GetParError(1);
744 data[key].t0fittype = 2;
745 return tresHist->GetFunction("gaus")->GetParameter(1);
746 }
747
748 tresHist->Fit("gaus","QR","",tresHist->GetFunction("gaus")->GetParameter(1) - tresHist->GetFunction("gaus")->GetParameter(2),tresHist->GetFunction("gaus")->GetParameter(1) + tresHist->GetFunction("gaus")->GetParameter(2));
749 if (tresHist->GetFunction("gaus")->GetParameter(1) - tresHist->GetFunction("gaus")->GetParameter(2) < m_mintres ||
750 tresHist->GetFunction("gaus")->GetParameter(1) + tresHist->GetFunction("gaus")->GetParameter(2) > m_maxtres) {
751 data[key].tresMean = mean;
752 data[key].t0err = rms;
753 data[key].t0fittype = 1;
754 return mean;
755 }
756
757 data[key].tres = tresHist->GetFunction("gaus")->GetParameter(2);
758 data[key].tresMean = tresHist->GetFunction("gaus")->GetParameter(1);
759 data[key].t0err = tresHist->GetFunction("gaus")->GetParError(1);
760 data[key].t0fittype = 2;
761
762 // protection against wrong fits:
763 if ( std::abs(tresHist->GetFunction("gaus")->GetParameter(2))>10 || std::abs(tresHist->GetFunction("gaus")->GetParameter(1))>10 ) {
764
765 data[key].tres = tresHist->GetRMS();
766 data[key].tresMean = tresHist->GetMean();
767 data[key].t0err = tresHist->GetRMS();
768 data[key].t0fittype = 6;
769 return tresHist->GetMean();
770 }
771
772
773 gStyle->SetOptFit(1111);
774
775 return tresHist->GetFunction("gaus")->GetParameter(1);
776 }
void mean(std::vector< double > &bins, std::vector< double > &values, const std::vector< std::string > &files, const std::string &histname, const std::string &tplotname, const std::string &label="")

◆ ATLAS_NOT_THREAD_SAFE() [3/3]

TDirectory *Calibrator::Calibrate ATLAS_NOT_THREAD_SAFE ( TDirectory * dir,
std::string key,
const std::string & opt,
caldata * caldata_above )

Definition at line 798 of file Calibrator.cxx.

798 { // Thread unsafe FitResidual, FitRt, FitTimeResidual are used.
799
800 //set some bool flags
801 bool calrt=opt.find('R')!=std::string::npos;
802 bool calt0=opt.find('T')!=std::string::npos;
803 bool donothing=opt.find('N')!=std::string::npos;
804 bool prnt=opt.find('P')!=std::string::npos;
805 bool useref=opt.find('B')!=std::string::npos;
806 bool useshortstw=opt.find('S')!=std::string::npos;
807
808 if (donothing) return dir;
809
810 data[key].calflag=true;
811
812 bool enough_t0 = data[key].ntres>=m_mint0stat;
813 bool enough_rt = data[key].nrt>=m_minrtstat;
814
815 if ((enough_rt && calrt) || (enough_t0 && calt0)) {
816 m_hdirs[key] = dir->mkdir(Form("%s%s",m_name.data(),key.data()));
817 m_hdirs[key]->cd();
818 if(level<2) {
819 std::cout << " Calibrator::Calibrate Parent Directory: " << std::string(dir->GetPath()) << std::endl;
820 std::cout << " Calibrator::Calibrate key " << key << " name " << m_name << " entries for rt " << data[key].nrt << " min req " << m_minrtstat << std::endl;
821 std::cout << " Calibrator::Calibrate key " << key << " entries for t0 " << data[key].ntres << " min req " << m_mint0stat << std::endl;
822 std::cout << " Calibrator::Calibrate Directory: " << std::string(m_hdirs[key]->GetPath()) << std::endl;
823 }
824 }
825 else m_hdirs[key]=dir;
826
827 //Fit also the residual if an rt or t0 calibration is made
828 if ((enough_rt && calrt) || (enough_t0 && calt0)) {
829
830 m_resHists[key] = new TH1F("residual","residual",m_nbinsres,m_minres,m_maxres);
831 int ncont=0;
832 for (int i=0;i<100;i++) {
833 m_resHists[key]->SetBinContent(i+1,data[key].reshist[i]);
834 ncont+=data[key].reshist[i];
835 }
836 m_resHists[key]->SetEntries((int)data[key].nres);
837 if (bequiet) m_resHists[key]->SetDirectory(nullptr);
838 if(ncont==0) {
839 std::cout << " No bin content " << std::endl;
840 } else {
841 FitResidual(key,m_resHists[key]);
842 }
843
844 }
845
846
847 if (prnt) printf("Calibrator::calibrate %8s %-14s: \n",m_name.data(),key.data());
848
849 //Calibrate r-m_t
850 if (nort){
851 //use old data
852 //std::cout << " use old data " << std::endl;
853 data[key].rtflag=true;
854 data[key].rtt0=caldata_above->rtt0;
855 data[key].rtgraph=caldata_above->rtgraph;
856 for (int i=0;i<4;i++) data[key].rtpar[i]=data[key].oldrtpar[i];
857 if (prnt) printf("RT << %7i (%8.1e) %8.1e %8.1e %8.1e, %3.2f : ", data[key].nrt,data[key].rtpar[0],data[key].rtpar[1],data[key].rtpar[2],data[key].rtpar[3], data[key].rtt0);
858 }
859 else{
860 //do fit
861 if ((calrt && enough_rt) || level==0) {
862 data[key].rtflag=true;
863 m_rtHists[key] = new TH2F("rt-relation","rt-relation",m_nbinst,m_mint,m_maxt,m_nbinsr,m_minr,m_maxr); //make a root histogram
864 for (int i=0;i<m_nbinst;i++) {
865 for (int j=0;j<m_nbinsr;j++) {
866 m_rtHists[key]->SetBinContent(i+1,j+1,data[key].rthist[i+m_nbinst*j]);
867 }
868 }
869 m_rtHists[key]->SetEntries(data[key].nrt);
870 if(level<2) std::cout << " Calibrator::calibrate do rt fit for key " << key << " with number of entries " << m_rtHists[key]->GetEntries() << std::endl;
871 if (bequiet) {
872 m_rtHists[key]->SetDirectory(nullptr);
873 if(level<2) std::cout << " remove this level from calout.root file " << std::endl;
874 }
875 data[key].rtt0=FitRt(key,opt,m_rtHists[key],m_hdirs[key]); //do the fit
876 if (prnt) printf("RT %7i (%8.1e) %8.1e %8.1e %8.1e, %3.2f : ", data[key].nrt, data[key].rtpar[0], data[key].rtpar[1], data[key].rtpar[2], data[key]. rtpar[3], data[key].rtt0);
877 }
878 else{
879 //use data from level above
880 data[key].rtgraph=caldata_above->rtgraph;
881 data[key].rtt0=caldata_above->rtt0;
882 for (int i=0;i<4;i++) data[key].rtpar[i]=caldata_above->rtpar[i];
883 if (prnt) printf("RT /\\ %7i (%8.1e) %8.1e %8.1e %8.1e, %3.2f : ", data[key].nrt,data[key].rtpar[0],data[key].rtpar[1],data[key].rtpar[2],data[key].rtpar[3], data[key].rtt0);
884 }
885 }
886
887
888 //Calibrate t0
889 if (not0){
890 //use old data
891 //std::cout << " use old data for T0 " std::endl;
892 data[key].t0flag=true;
893 data[key].t0=data[key].oldt02 + data[key].rtt0;
894 data[key].t0err=0;
895 data[key].t0off=0;
896 data[key].t0fittype = 5;
897 if (prnt) printf("T0 << %7i %05.2f%+05.2f%+05.2f=%05.2f \n", data[key].ntres, data[key].oldt02, 0.0, data[key].rtt0, data[key].t0);
898 }
899 else{
900 if (useref && level==5){
901 //use chip reference values
902 //std::cout << " use chip reference for T0 " << std::endl;
903 data[key].t0flag=true;
904 data[key].t0=caldata_above->t0 + data[key].reft0 + data[key].rtt0;
905 data[key].t0err=caldata_above->t0err;
906 data[key].t0off=data[key].t0-caldata_above->t0;
907 data[key].t0fittype = 3;
908 if (prnt) printf("T0 ** %7i %05.2f%+05.2f%+05.2f=%05.2f \n", data[key].ntres, caldata_above->t0, data[key].reft0, data[key].rtt0, data[key].t0);
909 }
910 else {
911 //do fit
912 if ((calt0 && enough_t0) || level==0) {
913 data[key].t0flag=true;
914 m_tresHists[key] = new TH1F("timeresidual","time residual",m_nbinstres,m_mintres,m_maxtres); //make a root histogram
915
916 for (int i=0;i<100;i++) {
917 m_tresHists[key]->SetBinContent(i+1,data[key].m_treshist[i]);
918 }
919 m_tresHists[key]->SetEntries(data[key].ntres);
920 if(level<2) std::cout << " Calibrator::calibrate do t0 fit for key " << key << " Numb. entries " << m_rtHists[key]->GetEntries() << std::endl;
921 if (bequiet) {
922 if(level<2) std::cout << " Calibrator::calibrate: remove this level from calibout.root " << std::endl;
923 m_tresHists[key]->SetDirectory(nullptr);
924 }
925
926 data[key].t0=data[key].oldt02 + FitTimeResidual(key,m_tresHists[key]) + data[key].rtt0 + m_t0shift; //do the fit and modify t0
927 data[key].t0off=data[key].t0-caldata_above->t0; //calculate t0 offset from level above
928 if (data[key].t0<0) data[key].t0=0;
929
930 if (prnt) printf("T0 %7i %05.2f%+05.2f%+05.2f%+05.2f=%05.2f \n", data[key].ntres, data[key].oldt02, data[key].t0-data[key].oldt02-data[key].rtt0, data[key].rtt0, m_t0shift, data[key].t0);
931
932
933 }
934 //use data from level above
935 else {
936 //std::cout << " use data from the level above " << std::endl;
937 //TEMP FIX to not destroy right T0s
938 if (data[key].oldt02 + (caldata_above->t0 - caldata_above->oldt02) >0) data[key].t0=data[key].oldt02 + (caldata_above->t0 - caldata_above->oldt02);
939 else data[key].t0= 0;
940 //TEMP FIX to not destroy right T0s
941
942 //add the short straw correction here. In this way, the shift is only done when contants at STRAW level come from level above.
943 if ((level == 6 && useshortstw) && std::abs(data[key].det)<2 && (data[key].lay==0 && data[key].stl<9) ) data[key].t0=caldata_above->t0-0.75;
944 data[key].t0err=caldata_above->t0err;
945 data[key].t0off=data[key].t0-caldata_above->t0 + data[key].rtt0;
946 data[key].t0fittype = 4;
947 if (prnt) printf("T0 /\\ %7i %05.2f%+05.2f%+05.2f=%05.2f \n", data[key].ntres, caldata_above->t0, caldata_above->oldt02, data[key].oldt02 , data[key].t0);
948 }
949 }
950 }
951
952
953 if (prnt && !bequiet) std::cout << " H";
954
955 if (prnt) std::cout << std::endl;
956
957 return m_hdirs[key];
958
959}
float reshist[100]
the 1D residual histogram (100 bins)
Definition Calibrator.h:157
float rthist[640]
the 2D rt histogram (20x32 bins)
Definition Calibrator.h:158
static Double_t t0
float t0
the new t0
Definition Calibrator.h:126
float t0err
the new to error
Definition Calibrator.h:127
float rtt0
the t0 fron the R-t fit
Definition Calibrator.h:130
RtGraph * rtgraph
the rt graph
Definition Calibrator.h:148
std::array< float, 4 > rtpar
the new rt-parameters
Definition Calibrator.h:132
double oldt02
the old sub-module t0 (average of t0 for all straws in the module)
Definition Calibrator.h:137
TH2F(name, title, nxbins, bins_par2, bins_par3, bins_par4, bins_par5=None, bins_par6=None, path='', **kwargs)
TH1F(name, title, nxbins, bins_par2, bins_par3=None, path='', **kwargs)

◆ pol3deg()

double pol3deg ( double * x,
double * par )

Definition at line 303 of file Calibrator.cxx.

303 {
304 double r = x[0];
305 double t = par[0]+r*(par[1]+(r*(par[2]+r*par[3])));
306 return t;
307}
#define x
int r
Definition globals.cxx:22

◆ pol3deg2()

double pol3deg2 ( double * x,
double * par )

Definition at line 309 of file Calibrator.cxx.

309 {
310 double t = x[0];
311 double r = par[1]*(t-par[0]) + par[2]*(t-par[0])*(t-par[0]) + par[3]*(t-par[0])*(t-par[0])*(t-par[0]);
312 return r;
313}

◆ rtrel_dines()

double rtrel_dines ( double * x,
double * par )

Definition at line 326 of file Calibrator.cxx.

326 {
327 double rmin0 = par[0];
328 double rho = par[1];
329 double v = par[2];
330 double t_const = par[3];
331 double t = x[0];
332 double r_squared = (4*rho*rho*std::sin(v*(t-t_const)/(2*rho))-rmin0*rmin0)/(1-0.25*rmin0*rmin0);
333 double r=r_squared>0 ? std::sqrt(r_squared) : 0.0;
334 return r;
335}

◆ trrel_dines()

double trrel_dines ( double * x,
double * par )

Definition at line 316 of file Calibrator.cxx.

316 {
317 double rmin0 = par[0];
318 double rho = par[1];
319 double v = par[2];
320 double t_const = par[3];
321 double r = x[0];
322 double t = t_const + 2*rho/v*std::asin(std::sqrt(rmin0*rmin0*(1-0.25*r*r)+r*r)/(2*rho));
323 return t;
324}