Skip to content

Commit

Permalink
Merge pull request #257 from VERITAS-Observatory/v490.6-fitstat-cuts
Browse files Browse the repository at this point in the history 
v490.6 - add additional cuts in mscw for event quality
  • Loading branch information
@GernotMaier
GernotMaier authored Jan 31, 2024
2 parents 32f178f + afb62eb commit d572b38
Show file tree
Hide file tree
Showing 5 changed files with 65 additions and 24 deletions.
14 changes: 10 additions & 4 deletions inc/VDispAnalyzer.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,8 @@ class VDispAnalyzer
float fdistance_max;
float floss_max;
float fFui_min;
float fWidth_min;
float fFitstat_min;
bool fDispErrorWeighting;
float fDispErrorExponential;

Expand Down Expand Up @@ -109,7 +111,8 @@ class VDispAnalyzer
double* iR, double iEHeight,
double iMCEnergy,
double* img_fui,
float* img_pedvar );
float* img_pedvar,
int* img_fitstat );

void calculateMeanDirection( float& xs, float& ys,
vector< float >& x, vector< float >& y,
Expand All @@ -134,7 +137,8 @@ class VDispAnalyzer
double* img_fui,
float* img_pedvar,
double* pointing_dx, double* pointing_dy,
bool UseIntersectForHeadTail );
bool UseIntersectForHeadTail,
int* img_fitstat );

vector< float > calculateExpectedDirectionError_or_Sign( unsigned int i_ntel, float iArrayElevation, float iArrayAzimuth,
ULong64_t* iTelType,
Expand All @@ -144,7 +148,7 @@ class VDispAnalyzer
double* img_tgrad, double* img_loss, int* img_ntubes,
double* img_weight,
double xoff_4, double yoff_4,
double* img_fui, float* img_pedvar );
double* img_fui, float* img_pedvar, int* img_fitstat );

float evaluate( float iWidth, float iLength, float iAsymm, float iDist,
float iSize, float iPedvar, float itgrad, float iLoss,
Expand Down Expand Up @@ -241,13 +245,15 @@ class VDispAnalyzer
}
void setQualityCuts( unsigned int iNImages_min = 0, float iAxesAngles_min = 0.,
float imaxdist = 1.e5, float imaxloss = 1.,
float iminfui = 0. )
float iminfui = 0., float iminwidth = -1., float iminfitstat = -10. )
{
fAxesAngles_min = iAxesAngles_min;
fNImages_min = iNImages_min;
fdistance_max = imaxdist;
floss_max = imaxloss;
fFui_min = iminfui;
fWidth_min = iminwidth;
fFitstat_min = iminfitstat;
}
void setTelescopeTypeList( vector<ULong64_t> iTelescopeTypeList );
void setTelescopeFOV( vector< float > iTelFOV )
Expand Down
4 changes: 3 additions & 1 deletion inc/VTableLookupRunParameter.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -79,6 +79,8 @@ class VTableLookupRunParameter : public TNamed, public VGlobalRunParameter
double fmaxloss; // note: same for all telescope types
double fminfui ; // note: same for all telescope types
double fminsize; // note: same for all telescope times
double fminwidth; // note: same for all telescope times
double fminfitstat; // note: same for all telescope times
// seed for random selection of showers before table filling
double fSelectRandom;
int fSelectRandomSeed;
Expand Down Expand Up @@ -110,6 +112,6 @@ class VTableLookupRunParameter : public TNamed, public VGlobalRunParameter
void print( int iB = 0 );
void printHelp();

ClassDef( VTableLookupRunParameter, 30 );
ClassDef( VTableLookupRunParameter, 31 );
};
#endif
27 changes: 18 additions & 9 deletions src/VDispAnalyzer.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -462,7 +462,8 @@ void VDispAnalyzer::calculateMeanDispDirection( unsigned int i_ntel,
float* img_pedvar,
double* pointing_dx,
double* pointing_dy,
bool UseIntersectForHeadTail )
bool UseIntersectForHeadTail,
int* img_fitstat )
{
// reset values from previous event
f_disp = -99.;
Expand All @@ -477,7 +478,7 @@ void VDispAnalyzer::calculateMeanDispDirection( unsigned int i_ntel,
|| !img_asym || !img_tgrad
|| !img_loss || !img_ntubes
|| !img_weight || !img_fui
|| !img_pedvar )
|| !img_pedvar || !img_fitstat )
{
return;
}
Expand All @@ -501,7 +502,9 @@ void VDispAnalyzer::calculateMeanDispDirection( unsigned int i_ntel,
if( img_size[i] > 0. && img_length[i] > 0. && img_ntubes[i] > 0
&& sqrt( img_cen_x[i]*img_cen_x[i] + img_cen_y[i]*img_cen_y[i] ) < fdistance_max
&& img_loss[i] < floss_max
&& img_fui[i] > fFui_min )
&& img_fui[i] > fFui_min
&& img_width[i] > fWidth_min
&& ( img_fitstat[i] < 1 || img_fitstat[i] >= fFitstat_min ) )
{
disp = evaluate( ( float )img_width[i], ( float )img_length[i], ( float )img_asym[i],
( float )sqrt( img_cen_x[i] * img_cen_x[i] + img_cen_y[i] * img_cen_y[i] ),
Expand Down Expand Up @@ -597,7 +600,8 @@ vector< float > VDispAnalyzer::calculateExpectedDirectionError_or_Sign( unsigned
double xoff_4,
double yoff_4,
double* img_fui,
float* img_pedvar )
float* img_pedvar,
int* img_fitstat )
{
vector< float > i_disp( i_ntel, -9999. );
// make sure that all data arrays exist
Expand All @@ -607,7 +611,7 @@ vector< float > VDispAnalyzer::calculateExpectedDirectionError_or_Sign( unsigned
|| !img_asym || !img_tgrad
|| !img_loss || !img_ntubes
|| !img_weight || !img_fui
|| !img_pedvar )
|| !img_pedvar || !img_fitstat )
{
return i_disp;
}
Expand All @@ -621,7 +625,9 @@ vector< float > VDispAnalyzer::calculateExpectedDirectionError_or_Sign( unsigned
if( img_size[i] > 0. && img_length[i] > 0.
&& sqrt( img_cen_x[i]*img_cen_x[i] + img_cen_y[i]*img_cen_y[i] ) < fdistance_max
&& img_loss[i] < floss_max
&& img_fui[i] > fFui_min )
&& img_fui[i] > fFui_min
&& img_width[i] > fWidth_min
&& ( img_fitstat[i] < 1 || img_fitstat[i] >= fFitstat_min ) )
{
i_disp[i] = evaluate( ( float )img_width[i], ( float )img_length[i], ( float )img_asym[i],
( float )sqrt( img_cen_x[i] * img_cen_x[i] + img_cen_y[i] * img_cen_y[i] ),
Expand Down Expand Up @@ -706,7 +712,8 @@ void VDispAnalyzer::calculateEnergies( unsigned int i_ntel,
double iEHeight,
double iMCEnergy,
double* img_fui,
float* img_pedvar )
float* img_pedvar,
int* img_fitstat )
{
fdisp_energy = -99.;
fdisp_energy_chi = 1.e-10;
Expand All @@ -723,7 +730,7 @@ void VDispAnalyzer::calculateEnergies( unsigned int i_ntel,
|| !img_asym || !img_tgrad
|| !img_loss || !img_ntubes
|| !img_weight || !iRcore
|| !img_fui )
|| !img_fui || !img_fitstat )
{
return;
}
Expand All @@ -738,7 +745,9 @@ void VDispAnalyzer::calculateEnergies( unsigned int i_ntel,
if( img_size[i] > 0. && iRcore[i] > 0. && iArrayElevation > 0.
&& sqrt( img_cen_x[i]*img_cen_x[i] + img_cen_y[i]*img_cen_y[i] ) < fdistance_max
&& img_loss[i] < floss_max
&& img_fui[i] > fFui_min )
&& img_fui[i] > fFui_min
&& img_width[i] > fWidth_min
&& ( img_fitstat[i] < 1 || img_fitstat[i] >= fFitstat_min ) )
{
fdisp_energy_T[i] = fTMVADispAnalyzer->evaluate(
( float )img_width[i], ( float )img_length[i],
Expand Down
24 changes: 15 additions & 9 deletions src/VTableLookupDataHandler.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -561,6 +561,7 @@ int VTableLookupDataHandler::fillNextEvent( bool bShort )
}
}
fCurrentNoiseLevel[i] = ftpars[i]->meanPedvar_Image;
fFitstat[i] = ftpars[i]->Fitstat;
if( !bShort )
{
fmeanPedvar_ImageT[i] = ftpars[i]->meanPedvar_Image;
Expand All @@ -577,7 +578,6 @@ int VTableLookupDataHandler::fillNextEvent( bool bShort )
fmaxindex2[i] = ftpars[i]->index_of_max[1];
fmaxindex3[i] = ftpars[i]->index_of_max[2];
ftchisq_x[i] = ftpars[i]->tchisq_x;
fFitstat[i] = ftpars[i]->Fitstat;
}
}
else
Expand Down Expand Up @@ -616,7 +616,9 @@ int VTableLookupDataHandler::fillNextEvent( bool bShort )
fDispAnalyzerEnergy->setQualityCuts( fSSR_NImages_min, fSSR_AxesAngles_min,
fTLRunParameter->fmaxdist,
fTLRunParameter->fmaxloss,
fTLRunParameter->fminfui );
fTLRunParameter->fminfui,
fTLRunParameter->fminwidth,
fTLRunParameter->fminfitstat );
fDispAnalyzerEnergy->calculateEnergies(
getNTel(),
fArrayPointing_Elevation, fArrayPointing_Azimuth,
Expand All @@ -632,7 +634,7 @@ int VTableLookupDataHandler::fillNextEvent( bool bShort )
getDistanceToCoreTel(),
fEmissionHeightMean,
fMCEnergy,
ffui, fmeanPedvar_ImageT );
ffui, fmeanPedvar_ImageT, fFitstat );

// fill results
setEnergy( fDispAnalyzerEnergy->getEnergy(), false );
Expand Down Expand Up @@ -718,7 +720,7 @@ void VTableLookupDataHandler::doStereoReconstruction()
floss, fntubes,
getWeight(),
fXoff_intersect, fYoff_intersect,
ffui, fmeanPedvar_ImageT );
ffui, fmeanPedvar_ImageT, fFitstat );
}
////////////////////////////////////////////////////////////////////
// estimate disp head/tail sign
Expand All @@ -738,7 +740,7 @@ void VTableLookupDataHandler::doStereoReconstruction()
floss, fntubes,
getWeight(),
fXoff_intersect, fYoff_intersect,
ffui, fmeanPedvar_ImageT );
ffui, fmeanPedvar_ImageT, fFitstat );
}

// use weighting calculated from disp error
Expand All @@ -747,7 +749,9 @@ void VTableLookupDataHandler::doStereoReconstruction()
fDispAnalyzerDirection->setQualityCuts( fSSR_NImages_min, fSSR_AxesAngles_min,
fTLRunParameter->fmaxdist,
fTLRunParameter->fmaxloss,
fTLRunParameter->fmaxloss );
fTLRunParameter->fminfui,
fTLRunParameter->fminwidth,
fTLRunParameter->fminfitstat );
fDispAnalyzerDirection->setTelescopeFOV( fTelFOV );
fDispAnalyzerDirection->calculateMeanDispDirection(
getNTel(),
Expand All @@ -764,7 +768,7 @@ void VTableLookupDataHandler::doStereoReconstruction()
iDispError, iDispSign,
ffui, fmeanPedvar_ImageT,
fpointing_dx, fpointing_dy,
fTLRunParameter->fDisp_UseIntersectForHeadTail );
fTLRunParameter->fDisp_UseIntersectForHeadTail, fFitstat );
// reconstructed direction by disp method:
fXoff = fDispAnalyzerDirection->getXcoordinate_disp();
fYoff = fDispAnalyzerDirection->getYcoordinate_disp();
Expand Down Expand Up @@ -1487,6 +1491,8 @@ bool VTableLookupDataHandler::setOutputFile( string iOutput, string iOption, str
fOTree->Branch( "ntubes", fntubes, iTT );
sprintf( iTT, "nsat[%d]/s", fNTel );
fOTree->Branch( "nsat", fnsat, iTT );
sprintf( iTT, "fui[%d]/D", fNTel );
fOTree->Branch( "fui", ffui, iTT );
sprintf( iTT, "nlowgain[%d]/s", fNTel );
fOTree->Branch( "nlowgain", fnlowgain, iTT );
sprintf( iTT, "alpha[%d]/D", fNTel );
Expand All @@ -1507,9 +1513,9 @@ bool VTableLookupDataHandler::setOutputFile( string iOutput, string iOption, str
fOTree->Branch( "tgrad_x", ftgrad_x, iTT );
sprintf( iTT, "tchisq_x[%d]/D", fNTel );
fOTree->Branch( "tchisq_x", ftchisq_x, iTT );
sprintf( iTT, "Fitstat[%d]/I", fNTel );
fOTree->Branch( "Fitstat", fFitstat, iTT );
}
sprintf( iTT, "Fitstat[%d]/I", fNTel );
fOTree->Branch( "Fitstat", fFitstat, iTT );
fOTree->Branch( "DispNImages", &fnxyoff, "DispNImages/i" );
fOTree->Branch( "DispXoff_T", fXoff_T, "DispXoff_T[NImages]/F" );
fOTree->Branch( "DispYoff_T", fYoff_T, "DispYoff_T[NImages]/F" );
Expand Down
20 changes: 19 additions & 1 deletion src/VTableLookupRunParameter.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -38,6 +38,8 @@ VTableLookupRunParameter::VTableLookupRunParameter()
fmaxdist = 50000.;
fmaxloss = 1.;
fminfui = 0.;
fminwidth = -1.;
fminfitstat = -10.;
fSelectRandom = -1.;
fSelectRandomSeed = 17;
fRerunStereoReconstruction = false;
Expand Down Expand Up @@ -294,6 +296,14 @@ bool VTableLookupRunParameter::fillParameters( int argc, char* argv[] )
{
fminfui = atof( iTemp.substr( iTemp.rfind( "=" ) + 1, iTemp.size() ).c_str() );
}
else if( iTemp.find( "-minwidth" ) < iTemp.size() )
{
fminwidth = atof( iTemp.substr( iTemp.rfind( "=" ) + 1, iTemp.size() ).c_str() );
}
else if( iTemp.find( "-minfitstat" ) < iTemp.size() )
{
fminfitstat = atof( iTemp.substr( iTemp.rfind( "=" ) + 1, iTemp.size() ).c_str() );
}
else if( iTemp.find( "-maxdistancetocameracenter" ) < iTemp.size() )
{
fMC_distance_to_cameracenter_max = atof( iTemp.substr( iTemp.rfind( "=" ) + 1, iTemp.size() ).c_str() );
Expand Down Expand Up @@ -564,7 +574,15 @@ void VTableLookupRunParameter::print( int iP )
}
if( fminfui > 0. )
{
cout << "\t BDT TMVA stereo reconstruction fui cut < " << fminfui << endl;
cout << "\t BDT TMVA stereo reconstruction fui cut > " << fminfui << endl;
}
if( fminwidth > 0. )
{
cout << "\t BDT TMVA stereo reconstruction min width cut < " << fminwidth << endl;
}
if( fminfitstat > 0. )
{
cout << "\t BDT TMVA stereo reconstruction min fitstat cut < " << fminfitstat << endl;
}
cout << "\t Head/tail uncertainty: ";
if( fDisp_UseIntersectForHeadTail )
Expand Down

0 comments on commit d572b38

Please sign in to comment.