00001 #include "StEEmcPointUtils.h"
00002 #include <algorithm>
00003
00004
00005 StEEmcPointUtils::StEEmcPointUtils()
00006 {
00007 }
00008
00009
00010 Bool_t StEEmcPointUtils::opening_angle_cut ( StEEmcPoint &p1, StEEmcPoint &p2, Float_t energy, Float_t mass, Float_t frac )
00011 {
00012 Float_t oangle = p1.position().Angle( p2.position() );
00013 Float_t oangle_min = 2.0 * mass * frac / energy;
00014 return (oangle>oangle_min);
00015 }
00016
00017
00018 StEEmcPointVec_t StEEmcPointUtils::find_best_points( StEEmcPointVec_t points, Float_t em_energy, Bool_t &stat, Float_t &chi2, Int_t npoint )
00019 {
00020
00023 chi2=1.0E9;
00024 if ( npoint > (Int_t)points.size() ) {
00025 stat=false;
00026 return points;
00027 }
00028
00031 Bool_t go=true;
00032 Float_t chi2_min=1.0E9;
00033 StEEmcPointVec_t good_perm;
00034 while ( go ) {
00035
00038 if ( no_parallel_smd_clusters( points, npoint ) )
00039 for ( Int_t ipoint=0; ipoint < npoint; ipoint++ ) {
00040 Float_t edep_predict = 2.0 * kEEmcSmdSF * em_energy;
00041 Float_t ediff1 = points[ipoint].energy() - edep_predict;
00042 Float_t ediff2 = points[ipoint].cluster(0).energy() - points[ipoint].cluster(1).energy();
00043 ediff2 *= 2.0;
00044 Float_t mychi2=ediff1*ediff1 + ediff2*ediff2;
00045 if ( mychi2 < chi2_min ) {
00046 good_perm=points;
00047 chi2_min=mychi2;
00048 }
00049
00050 }
00051
00052 go = std::next_permutation ( points.begin(), points.end() );
00053
00054 }
00055
00056 return points;
00057
00058 }
00059
00060
00061 Bool_t StEEmcPointUtils::no_parallel_smd_clusters( StEEmcPointVec_t points, Int_t npoint )
00062 {
00063
00064 if ( npoint <= 1 ) return true;
00065
00066 for ( Int_t i = 0; i < npoint-1; i++ ) {
00067 Float_t u1=points[i].cluster(0).mean();
00068 Float_t v1=points[i].cluster(1).mean();
00069 for ( Int_t j = i+1; j < npoint; j++ ) {
00070 Float_t u2=points[j].cluster(0).mean();
00071 Float_t v2=points[j].cluster(1).mean();
00072 if ( u1==u2 ) return false;
00073 if ( v1==v2 ) return false;
00074 }
00075 }
00076
00078 return true;
00079
00080 }
