Public Member Functions | Public Attributes | Private Attributes | List of all members
fuzz::FuzzyKMeanAlg Class Reference

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N21-01-15/FuzzyKVertex/FuzzyKMeanAlg.h"

Public Member Functions

 FuzzyKMeanAlg (const FuzzyKMeanParams &params)
 
 ~FuzzyKMeanAlg ()
 
void LoadCluster (rb::Cluster slice)
 Function to take in a cluster and build the hitlist needed to do the clustering. More...
 
void DoClustering (unsigned int minClust, unsigned int maxClust, double xy, double z, geo::View_t view)
 A function to control all steps in clustering in 2D on a group of hits. More...
 
void MakeClusters (unsigned int n)
 A function to calculate membership weights for a specified number of clusters. More...
 
void SeedA (unsigned int n)
 A function to calculate the seeds used for the first round of clustering. More...
 
bool FindUnclustered ()
 A function to collect the hits that have not made it to a cluster yet and assign new cluster seeds based on that. More...
 
double Distance (double mu, double x, double dx)
 Calculate and return the distance from a hit to a cluster center. More...
 
void ComputeDij ()
 Compute all distances for j hits to i cluster centers. More...
 
void ComputeUij ()
 Compute the membership values for j hits to i cluster centers. More...
 
bool ComputeA ()
 Update the mean cluster center positions. More...
 
double Jpfuz ()
 Calculate and return the score evaluating the clustering. More...
 
bool CheckUnique ()
 Check if each of the cluster centers (after the algorith has converged) is unique, remove duplicates to reseed. More...
 
void MakeAngles (double vXY, double vZ, geo::View_t view)
 Function to turn the 2D hit coordinates into angles (radians) with respect to the vertex for clustering. More...
 
void MergeClusters ()
 Function to merge cluster centers together that have significant overlap. More...
 

Public Attributes

rb::HitList fHL
 List of hits in slice. More...
 
std::vector< double > fDist
 Distance from vertex. More...
 
std::vector< double > fW
 Hit weight, currently PE. More...
 
std::vector< double > fX
 Angle to vertex. More...
 
std::vector< double > fdX
 Uncertainty in angle. More...
 
std::vector< std::vector< double > > fD
 Dist clust i to sub. j. More...
 
std::vector< std::vector< double > > fU
 Assign cluster i to sub. j. More...
 
std::vector< double > fA
 Cluster centers. More...
 
std::vector< double > fSeed
 List of used seeds. More...
 
std::vector< std::pair< double, double > > fUnused
 Unused cell pos., uncert. More...
 
std::vector< std::pair< double, double > > fLMax
 List of max densities. More...
 

Private Attributes

FuzzyKMeanParams fParams
 

Detailed Description

Definition at line 55 of file FuzzyKMeanAlg.h.

Constructor & Destructor Documentation

fuzz::FuzzyKMeanAlg::FuzzyKMeanAlg ( const FuzzyKMeanParams params)
explicit

Definition at line 18 of file FuzzyKMeanAlg.cxx.

19  : fParams(params)
20  {
21  }
FuzzyKMeanParams fParams
fuzz::FuzzyKMeanAlg::~FuzzyKMeanAlg ( )

Definition at line 25 of file FuzzyKMeanAlg.cxx.

26  {
27  }

Member Function Documentation

bool fuzz::FuzzyKMeanAlg::CheckUnique ( )

Check if each of the cluster centers (after the algorith has converged) is unique, remove duplicates to reseed.

Function to check the final cluster centers and see if any are duplicates If a duplicate is found a replacement seed will be substituted or the duplicate seed will be removed if no replacement is available. Returns true if the clustering algorithm needs to run again with the updated seeds. Returns false if all clusters found are unique.

Definition at line 430 of file FuzzyKMeanAlg.cxx.

References fuzz::FuzzyKMeanParams::Duplicate, fA, stan::math::fabs(), fLMax, fParams, fSeed, MECModelEnuComparisons::i, calib::j, m, and moon_position_table_new3::second.

Referenced by DoClustering().

431  {
432 
433  bool recluster = false;
434  for (unsigned int i=0; i<fA.size(); i++){
435  for (unsigned int j=i+1; j<fA.size(); j++){
436  //check cluster centers for duplicates
437  if (fabs(fA[i]-fA[j]) < fParams.Duplicate()){
438  recluster = true;
439 
440  //if there is a duplicate see if a new seed on the list can be found
441  bool replace = true;
442  for (unsigned int k=0; k<fLMax.size(); k++){
443  replace = true;
444  for(unsigned int m=0; m<fSeed.size(); m++){
445  if (fSeed[m] == fLMax[k].second){
446  replace = false;
447  break;
448  }
449  }
450  if (replace){
451  fA[j] = fLMax[k].second;
452  fSeed.push_back(fLMax[k].second);
453  break;
454  }
455  }
456 
457  //if a replacement seed could not be found, remove the duplicate cluster
458  if (!replace){
459  fA.erase(fA.begin()+j);
460  j--;
461  }
462 
463  }//end if statement
464  }// end loop over j
465  }//end loop over i
466  return recluster;
467  }
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
std::vector< double > fSeed
List of used seeds.
std::vector< std::pair< double, double > > fLMax
List of max densities.
std::vector< double > fA
Cluster centers.
const double j
Definition: BetheBloch.cxx:29
FuzzyKMeanParams fParams
Atom< double > Duplicate
Definition: FuzzyKMeanAlg.h:50
bool fuzz::FuzzyKMeanAlg::ComputeA ( )

Update the mean cluster center positions.

Function to update cluster centers. Returns true if all the centers have stabalized within tolerance, false if another iteration is needed.

Definition at line 279 of file FuzzyKMeanAlg.cxx.

References om::cout, allTimeWatchdog::endl, fA, stan::math::fabs(), fdX, fParams, fU, fuzz::FuzzyKMeanParams::Fuzz, fX, MECModelEnuComparisons::i, calib::j, M_PI, cet::pow(), fuzz::FuzzyKMeanParams::PrintLevel, fuzz::FuzzyKMeanParams::Tol, and w.

Referenced by MakeClusters().

280  {
281  unsigned int i, j;
282  double s;
283  double w;
284  double musave;
285  bool done = true;
286  for (i=0; i<fA.size(); ++i) {
287  musave = fA[i];
288  s = 0.0;
289  w = 0.0;
290 
291  for (j=0; j<fX.size(); ++j) {
292  double dtheta = fX[j]-musave;
293  while (dtheta<-M_PI) dtheta += 2*M_PI;
294  while (dtheta> M_PI) dtheta -= 2*M_PI;
295  s += pow(fU[i][j],fParams.Fuzz())*dtheta/(fdX[j]*fdX[j]);
296  w += pow(fU[i][j],fParams.Fuzz())/(fdX[j]*fdX[j]);
297  }
298 
299  double thetaNew = musave + s/w;
300  while (thetaNew<-M_PI) thetaNew += 2*M_PI;
301  while (thetaNew> M_PI) thetaNew -= 2*M_PI;
302  fA[i] = thetaNew;
303  if (fabs(fA[i]-musave) > fParams.Tol()) done = false;
304  if (fParams.PrintLevel() >= 2) std::cout << "COMP MU["<<i<<"]="<<fA[i]<<std::endl;
305  }
306  return done;
307  }
Atom< double > Tol
Definition: FuzzyKMeanAlg.h:33
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
Atom< double > Fuzz
Definition: FuzzyKMeanAlg.h:32
constexpr T pow(T x)
Definition: pow.h:75
std::vector< double > fdX
Uncertainty in angle.
#define M_PI
Definition: SbMath.h:34
std::vector< double > fX
Angle to vertex.
const XML_Char * s
Definition: expat.h:262
std::vector< double > fA
Cluster centers.
std::vector< std::vector< double > > fU
Assign cluster i to sub. j.
const double j
Definition: BetheBloch.cxx:29
FuzzyKMeanParams fParams
OStream cout
Definition: OStream.cxx:6
Float_t w
Definition: plot.C:20
void fuzz::FuzzyKMeanAlg::ComputeDij ( )

Compute all distances for j hits to i cluster centers.

Function to calculate distance metric for all i clusters and j angles, with uncertainty.

Definition at line 236 of file FuzzyKMeanAlg.cxx.

References om::cout, d, Distance(), allTimeWatchdog::endl, fA, fD, fdX, fParams, fX, MECModelEnuComparisons::i, calib::j, and fuzz::FuzzyKMeanParams::PrintLevel.

Referenced by MakeClusters().

237  {
238  unsigned int i, j;
239  double d;
240  for (i=0; i<fA.size(); ++i) {
241  for (j=0; j<fX.size(); ++j) {
242  d = this->Distance(fA[i], fX[j], fdX[j]);
243  fD[i][j] = d*d;
244  if (fParams.PrintLevel() >= 2){
245  std::cout << "D[" << i << "," << j << "]=("
246  << fX[j] << "-" << fA[i] << ")/"
247  << fdX[j] << "=" << fD[i][j]
248  << std::endl;
249  }
250  }//end loop over j
251  }//end loop over it
252  }
std::vector< std::vector< double > > fD
Dist clust i to sub. j.
double Distance(double mu, double x, double dx)
Calculate and return the distance from a hit to a cluster center.
std::vector< double > fdX
Uncertainty in angle.
std::vector< double > fX
Angle to vertex.
std::vector< double > fA
Cluster centers.
Float_t d
Definition: plot.C:236
const double j
Definition: BetheBloch.cxx:29
FuzzyKMeanParams fParams
OStream cout
Definition: OStream.cxx:6
void fuzz::FuzzyKMeanAlg::ComputeUij ( )

Compute the membership values for j hits to i cluster centers.

Function to compute membership of cells to clusters.

Definition at line 257 of file FuzzyKMeanAlg.cxx.

References fuzz::FuzzyKMeanParams::Beta, om::cout, e, allTimeWatchdog::endl, stan::math::exp(), fA, fD, fParams, fU, fuzz::FuzzyKMeanParams::Fuzz, fX, MECModelEnuComparisons::i, calib::j, fuzz::FuzzyKMeanParams::PrintLevel, and std::sqrt().

Referenced by MakeClusters().

258  {
259  unsigned int i, j;
260  for (j=0; j<fX.size(); ++j) {
261  for (i=0; i<fA.size(); ++i) {
262  fU[i][j] = exp(-(fParams.Fuzz()*sqrt((double)fA.size())*fD[i][j])/fParams.Beta());
263  if (fU[i][j] != fU[i][j]) fU[i][j] = 1e-50; //clean up NaN problem
264  if (fU[i][j] < 1e-50) fU[i][j] = 1e-50; //set a membership floor
265  if (fParams.PrintLevel() >= 2){
266  std::cout << "U[" << i << "," << j << "]="
267  << fU[i][j] << " | ";
268  }
269  } //end loop over i
270  if (fParams.PrintLevel() >= 2) std::cout<<std::endl;
271  } //end loop over j
272  }
Atom< double > Beta
Definition: FuzzyKMeanAlg.h:34
std::vector< std::vector< double > > fD
Dist clust i to sub. j.
Atom< double > Fuzz
Definition: FuzzyKMeanAlg.h:32
T sqrt(T number)
Definition: d0nt_math.hpp:156
std::vector< double > fX
Angle to vertex.
std::vector< double > fA
Cluster centers.
std::vector< std::vector< double > > fU
Assign cluster i to sub. j.
fvar< T > exp(const fvar< T > &x)
Definition: exp.hpp:10
const double j
Definition: BetheBloch.cxx:29
FuzzyKMeanParams fParams
OStream cout
Definition: OStream.cxx:6
Float_t e
Definition: plot.C:35
double fuzz::FuzzyKMeanAlg::Distance ( double  mu,
double  x,
double  dx 
)

Calculate and return the distance from a hit to a cluster center.

Function to calculate distance from an angle to cluster center.

Parameters
muMean position of cluster center (in radians)
xPosition of a hit (in radians)
dxThe uncertainty in hit position (in radians)

Definition at line 225 of file FuzzyKMeanAlg.cxx.

References d, dx, and M_PI.

Referenced by ComputeDij().

226  {
227  double d = x-mu;
228  while (d<-M_PI) d += 2*M_PI;
229  while (d> M_PI) d -= 2*M_PI;
230  return d/dx;
231  }
#define M_PI
Definition: SbMath.h:34
double dx[NP][NC]
Float_t d
Definition: plot.C:236
void fuzz::FuzzyKMeanAlg::DoClustering ( unsigned int  minClust,
unsigned int  maxClust,
double  xy,
double  z,
geo::View_t  view 
)

A function to control all steps in clustering in 2D on a group of hits.

Parameters
minClustThe minimum number of cluster candidates to make
maxClustThe maximum number of cluster candidates to make
xyThe x or y coordinate (depending on view) of the vertex to cluster around
zThe z coordinate of the vertex to cluster around
viewThe detector view (either kX or kY) of the hits to use in clustering. This must be a 2d view

Function performs clustering using a fuzzy possibilistic clustering approach. Takes arguements for the minimum number of clusters to start with and the maximum number of clusters to allow. Unused hits will be left unclustered. Other arguments are a vertex position and geometry view

Definition at line 44 of file FuzzyKMeanAlg.cxx.

References CheckUnique(), om::cout, allTimeWatchdog::endl, fA, FindUnclustered(), fParams, fSeed, fU, fX, calib::j, M_PI, MakeAngles(), MakeClusters(), MergeClusters(), fuzz::FuzzyKMeanParams::PrintLevel, and SeedA().

Referenced by fuzz::FuzzyKVertex::AddProng(), and fuzz::FuzzyKVertex::MakeProngs().

46  {
47  this->MakeAngles(xy, z, view); //turn hits into angles around vertex
48  fA.resize(minClust);
49  fSeed.clear();
50  this->SeedA(minClust); //produce intitial cluster seeds
51 
52  //loop to increase number of clusters
53  while (fA.size() <=maxClust){
54  //optimize cluster centers and assign membership for set cluster number
55  this->MakeClusters(fA.size());
56 
57  //now check that all clusters centers produced are unique
58  //remove duplicates and recluster
59  if (this->CheckUnique() == true) continue;
60 
61  //check to see if there are cells without cluster membership
62  //stop clustering if no unassigned hits are found
63  if (this->FindUnclustered() == false) break;
64  }//end while
65 
66  //remove a cluster if it excedes the limit
67  if (fA.size() > maxClust) fA.pop_back();
68 
69  //now check to see if any clusters should be merged
70  this->MergeClusters();
71 
72  //print statements for debugging
73  if (fParams.PrintLevel() >=2){
74  std::cout<<"****************************************"<<std::endl;
75  std::cout<<"Made: "<<fA.size()<<" clusters!"<<std::endl;
76  for(unsigned int j=0;j<fA.size();j++){
77  std::cout<<"Cluster "<<j<<": "<<fA[j]*180.0/M_PI<<std::endl;
78  }
79  for (unsigned int j=0; j<fX.size(); ++j) {
80  for (unsigned int k=0; k<fA.size(); ++k) {
81  std::cout << "U[" << k << "," << j << "]="
82  << fU[k][j] << " | ";
83  }
85  }
86  }//end debugging
87 
88  }
std::vector< double > fSeed
List of used seeds.
bool FindUnclustered()
A function to collect the hits that have not made it to a cluster yet and assign new cluster seeds ba...
#define M_PI
Definition: SbMath.h:34
std::vector< double > fX
Angle to vertex.
void MakeAngles(double vXY, double vZ, geo::View_t view)
Function to turn the 2D hit coordinates into angles (radians) with respect to the vertex for clusteri...
void SeedA(unsigned int n)
A function to calculate the seeds used for the first round of clustering.
std::vector< double > fA
Cluster centers.
std::vector< std::vector< double > > fU
Assign cluster i to sub. j.
const double j
Definition: BetheBloch.cxx:29
FuzzyKMeanParams fParams
z
Definition: test.py:28
OStream cout
Definition: OStream.cxx:6
bool CheckUnique()
Check if each of the cluster centers (after the algorith has converged) is unique, remove duplicates to reseed.
void MergeClusters()
Function to merge cluster centers together that have significant overlap.
void MakeClusters(unsigned int n)
A function to calculate membership weights for a specified number of clusters.
bool fuzz::FuzzyKMeanAlg::FindUnclustered ( )

A function to collect the hits that have not made it to a cluster yet and assign new cluster seeds based on that.

Function to find cells left out of clusters, use them to pick new cluster seed Return true if a new seed candidate is found and clustering should continue, return false if all hits are in clusters, or there are no new seeds to try.

Definition at line 330 of file FuzzyKMeanAlg.cxx.

References fuzz::FuzzyKMeanParams::Bins, confusionMatrixTree::count, om::cout, release_diff::diff, allTimeWatchdog::endl, stan::math::exp(), fA, fdX, fLMax, fParams, fSeed, fU, fUnused, fX, MECModelEnuComparisons::i, calib::j, M_PI, make_pair(), fuzz::FuzzyKMeanParams::Memb, fuzz::FuzzyKMeanParams::PrintLevel, site_stats_from_log::reverse, moon_position_table_new3::second, fillBadChanDBTables::step, and w.

Referenced by DoClustering().

331  {
332  bool fContinue = false; //flag to continue clustering
333  fUnused.clear();
334  bool unwanted;
335  int count=0;
336  const double step = (double)(2.0*M_PI)/fParams.Bins();
337 
338  //loop over hits check to see if there are leftovers
339  for (unsigned int j=0; j<fX.size(); j++){
340  unwanted = true;
341  for (unsigned int i=0; i<fA.size(); i++){
342  //minimum membership set in FHiCL, break if a hit has membership somewhere
343  if (fU[i][j] >= fParams.Memb()){
344  unwanted = false;
345  break;
346  }
347  } //end loop over i
348  if (unwanted){
349  fUnused.push_back(std::make_pair(fX[j],fdX[j]));
350  count++;
351  }
352  } //end loop over j
353 
354  //if there are unclustered hits, find the most appropriate seed to add
355  if (count > 0){
356  std::vector<double> w(fParams.Bins());
357  double diff;
358  for (unsigned int k=0; k<w.size(); k++){
359  w[k] = 0.0;
360  for (unsigned int j=0; j<fUnused.size(); j++){
361  diff = -M_PI + k*step - fUnused[j].first;
362  while (diff < -M_PI) diff += 2*M_PI;
363  while (diff > M_PI) diff -= 2*M_PI;
364  w[k] += exp(-diff*diff/(fUnused[j].second*fUnused[j].second));
365  }
366  }
367 
368  //now find local maxima in weight distribution
369  fLMax.clear();
370  for (unsigned int k=0; k<w.size(); k++) {
371  if (k==0){
372  if ((w[k] >= w.back()) && (w[k] > w[k+1])) {
373  fLMax.push_back(std::make_pair(w[k],-M_PI+k*step));
374  }
375  }
376  else if (k==w.size()-1) {
377  if ((w[k] >= w[k-1]) && (w[k] > w.front())) {
378  fLMax.push_back(std::make_pair(w[k],-M_PI+k*step));
379  }
380  }
381  else {
382  if ((w[k] >= w[k-1]) && (w[k] > w[k+1])) {
383  fLMax.push_back(std::make_pair(w[k],-M_PI+k*step));
384  }
385  }
386  }
387  std::sort(fLMax.begin(),fLMax.end());
388  std::reverse(fLMax.begin(),fLMax.end());
389 
390  //loop over peaks in the unclustered distribution,
391  //try to add highest peak to the list of seed candidates
392  //In some cases an unclustered hit is at an angle where even
393  //if it has a seed the cluster will converge on an already found
394  //cluster center. So check that the seed to be added has not already
395  //been used
396  for (unsigned int i=0; i<fLMax.size(); i++){
397  bool addSeed = true;
398  //loop over seeds to see if it exists
399  for (unsigned int j=0; j<fSeed.size(); j++){
400  if (fSeed[j] == fLMax[i].second){
401  addSeed = false;
402  break;
403  }
404  }
405 
406  if (addSeed){
407  fA.push_back(fLMax[i].second);
408  fSeed.push_back(fLMax[i].second);
409  fContinue = true;
410  if (fParams.PrintLevel() >= 2){
411  std::cout<<"ADDING NEW CLUSTER AT: "
412  <<fLMax[i].second*180.0/M_PI<<std::endl;
413  }
414  break;
415  }
416 
417  }//end loop over i
418 
419  }//end if statement
420  return fContinue;
421  }
std::vector< double > fSeed
List of used seeds.
std::vector< double > fdX
Uncertainty in angle.
std::pair< Spectrum *, CheatDecomp * > make_pair(SpectrumLoaderBase &loader_data, SpectrumLoaderBase &loader_mc, HistAxis *axis, Cut *cut, const SystShifts &shift, const Var &wei)
Definition: DataMCLoad.C:336
#define M_PI
Definition: SbMath.h:34
std::vector< double > fX
Angle to vertex.
std::vector< std::pair< double, double > > fLMax
List of max densities.
std::vector< double > fA
Cluster centers.
std::vector< std::vector< double > > fU
Assign cluster i to sub. j.
fvar< T > exp(const fvar< T > &x)
Definition: exp.hpp:10
const double j
Definition: BetheBloch.cxx:29
FuzzyKMeanParams fParams
OStream cout
Definition: OStream.cxx:6
Atom< double > Memb
Definition: FuzzyKMeanAlg.h:35
std::vector< std::pair< double, double > > fUnused
Unused cell pos., uncert.
Float_t w
Definition: plot.C:20
double fuzz::FuzzyKMeanAlg::Jpfuz ( )

Calculate and return the score evaluating the clustering.

Function to calculate metric for clustering quality.

Definition at line 312 of file FuzzyKMeanAlg.cxx.

References fuzz::FuzzyKMeanParams::Beta, fA, fParams, fU, fuzz::FuzzyKMeanParams::Fuzz, fX, MECModelEnuComparisons::i, calib::j, cet::pow(), std::sqrt(), and sum.

313  {
314  unsigned int i, j;
315  double sum = 0.0;
316  for (i=0; i<fA.size(); ++i) {
317  for (j=0; j<fX.size(); ++j) {
318  sum += pow(fU[i][j],fParams.Fuzz());
319  }
320  }
321  sum = -sum*fParams.Beta()/(fParams.Fuzz()*fParams.Fuzz()*sqrt((double)fA.size()));
322  return sum;
323  }
Atom< double > Beta
Definition: FuzzyKMeanAlg.h:34
Atom< double > Fuzz
Definition: FuzzyKMeanAlg.h:32
T sqrt(T number)
Definition: d0nt_math.hpp:156
constexpr T pow(T x)
Definition: pow.h:75
std::vector< double > fX
Angle to vertex.
std::vector< double > fA
Cluster centers.
std::vector< std::vector< double > > fU
Assign cluster i to sub. j.
const double j
Definition: BetheBloch.cxx:29
FuzzyKMeanParams fParams
Double_t sum
Definition: plot.C:31
void fuzz::FuzzyKMeanAlg::LoadCluster ( rb::Cluster  slice)

Function to take in a cluster and build the hitlist needed to do the clustering.

Function to load in a cluster and build the hitlist.

Definition at line 32 of file FuzzyKMeanAlg.cxx.

References fHL, fParams, and fuzz::FuzzyKMeanParams::WeightByQ.

Referenced by fuzz::FuzzyKVertex::AddProng(), and fuzz::FuzzyKVertex::MakeProngs().

33  {
34  fHL.clear();
35  fHL = rb::HitList(slice, fParams.WeightByQ()); //initialize hitlist
36  }
std::vector< DAQHit > HitList
Definition: HitList.h:15
Atom< bool > WeightByQ
Definition: FuzzyKMeanAlg.h:31
rb::HitList fHL
List of hits in slice.
FuzzyKMeanParams fParams
void fuzz::FuzzyKMeanAlg::MakeAngles ( double  vXY,
double  vZ,
geo::View_t  view 
)

Function to turn the 2D hit coordinates into angles (radians) with respect to the vertex for clustering.

Parameters
vXYThe x or y (depending on view) vertex coordinate
vZThe z vertex coordinate
viewThe view, either kX or kY, of the hits to use. This needs to match the view supplied in vXY

Function to turn the hits in the hitlist from points into angles (in radians) with respect to a vertex candidate. This is done for a specified view.

Definition at line 129 of file FuzzyKMeanAlg.cxx.

References std::atan2(), fuzz::FuzzyKMeanParams::Cutoff, dist, fDist, fdX, fHL, fParams, fW, fX, MECModelEnuComparisons::i, fuzz::FuzzyKMeanParams::Par1, fuzz::FuzzyKMeanParams::Par2, fuzz::FuzzyKMeanParams::Par3, and util::pythag().

Referenced by DoClustering().

130  {
131  fX.clear();
132  fdX.clear();
133  fW.clear();
134  fDist.clear();
135 
136  //loop over hitlist calculating angles
137  for(unsigned int i=0;i<fHL.size();i++){
138  //make sure hit is in the right view
139  if (fHL[i].fView == view){
140  //calculate angle
141  fX.push_back(atan2(fHL[i].fXY-vXY,fHL[i].fZ-vZ));
142  //add weight (cell energy)
143  fW.push_back(fHL[i].fW);
144  //distance from vertex
145  double dist = util::pythag(fHL[i].fXY-vXY,fHL[i].fZ-vZ);
146  //if cell center is the same as the vertex, alter distance slightly to not divide by 0
147  if (dist == 0.0) dist = 0.0001;
148  fDist.push_back(dist);
149 
150  //uncertainty in angle, fit came from muon simulation
151  //could be adjusted for different degrees of fuzziness
152  //fit is of the form dx = [p1]/x + [p2] + [p3]*x
153  if (dist<fParams.Cutoff()) fdX.push_back(fParams.Par1()/dist + fParams.Par2() + fParams.Par3()*dist);
154  //beyond the cutoff distance keep uncertainty fixed
155  else fdX.push_back(fParams.Par1()/dist + fParams.Par2() + fParams.Par3()*fParams.Cutoff());
156  }//end if
157  }//end for loop
158 
159  }
std::vector< double > fDist
Distance from vertex.
std::vector< double > fdX
Uncertainty in angle.
Atom< double > Cutoff
Definition: FuzzyKMeanAlg.h:49
std::vector< double > fX
Angle to vertex.
double dist
Definition: runWimpSim.h:113
rb::HitList fHL
List of hits in slice.
FuzzyKMeanParams fParams
Atom< double > Par3
Definition: FuzzyKMeanAlg.h:42
std::vector< double > fW
Hit weight, currently PE.
double pythag(double x, double y)
2D Euclidean distance
Definition: MathUtil.h:29
Atom< double > Par2
Definition: FuzzyKMeanAlg.h:40
Atom< double > Par1
Definition: FuzzyKMeanAlg.h:38
T atan2(T number)
Definition: d0nt_math.hpp:72
void fuzz::FuzzyKMeanAlg::MakeClusters ( unsigned int  n)

A function to calculate membership weights for a specified number of clusters.

Parameters
nThe number of cluster centers to use

Function to make n clusters. Function loops until optimum number of cluster centers are found

Definition at line 94 of file FuzzyKMeanAlg.cxx.

References ComputeA(), ComputeDij(), ComputeUij(), om::cout, allTimeWatchdog::endl, fA, fD, fParams, fU, fX, MECModelEnuComparisons::i, M_PI, fuzz::FuzzyKMeanParams::MaxItr, getGoodRuns4SAM::n, and fuzz::FuzzyKMeanParams::PrintLevel.

Referenced by DoClustering().

95  {
96  unsigned int i;
97 
98  fD.resize(n); //resize vector holding distance to cluster centers
99  fU.resize(n); //resize vector holding membership probabilities
100  for (i=0; i<n; ++i) {
101  fD[i].resize(fX.size());
102  fU[i].resize(fX.size());
103  }
104 
105  if (fParams.PrintLevel() >= 2){
106  std::cout<<"Doing Clustering for "<<n<<" seeds: "<<std::endl;
107  for(unsigned int i=0; i<n; i++){
108  std::cout<<"seed["<<i<<"]: "<<fA[i]*180.0/M_PI<<std::endl;
109  }
110  }
111 
112  this->ComputeDij();//compute distances to cluster centers
113  this->ComputeUij();//compute membership probabilities
114 
115  //iterate cluster optimization until cluster centers stop changing
116  for (i=0; i<fParams.MaxItr(); ++i) {
117  bool done = this->ComputeA(); //calculate new clusters centers
118  this->ComputeDij(); //compute distances
119  this->ComputeUij(); //compute membership
120  if (done) break;
121  }
122  }
std::vector< std::vector< double > > fD
Dist clust i to sub. j.
bool ComputeA()
Update the mean cluster center positions.
#define M_PI
Definition: SbMath.h:34
std::vector< double > fX
Angle to vertex.
Atom< unsigned int > MaxItr
Definition: FuzzyKMeanAlg.h:36
void ComputeDij()
Compute all distances for j hits to i cluster centers.
std::vector< double > fA
Cluster centers.
std::vector< std::vector< double > > fU
Assign cluster i to sub. j.
FuzzyKMeanParams fParams
OStream cout
Definition: OStream.cxx:6
void ComputeUij()
Compute the membership values for j hits to i cluster centers.
void fuzz::FuzzyKMeanAlg::MergeClusters ( )

Function to merge cluster centers together that have significant overlap.

Function to check clusters and determine if they share enough hits to be merged into one cluster.

Definition at line 472 of file FuzzyKMeanAlg.cxx.

References release_diff::diff, fA, stan::math::fabs(), fD, fdX, fParams, util::frac(), fU, fuzz::FuzzyKMeanParams::Fuzz, fX, MECModelEnuComparisons::i, calib::j, M_PI, fuzz::FuzzyKMeanParams::Memb, fuzz::FuzzyKMeanParams::MergeAng, fuzz::FuzzyKMeanParams::MergeMemb, cet::pow(), gen_hdf5record::size, and w.

Referenced by DoClustering().

473  {
474  double diff;
475  unsigned int i=0;
476  unsigned int size = fA.size();
477 
478  while(i<size){
479  //loop over clusters and search for pairs within FHiCL specified cutoff (default of 30 degrees)
480  for (unsigned int j=i+1; j<fA.size(); j++){
481  diff = fabs(fA[i] - fA[j]);
482  while (diff > M_PI) diff -= 2*M_PI;
483 
484  if (diff < fParams.MergeAng()){
485  //now check out overlap between two clusters
486  std::vector<int> iClustCells; //cells with membership in i cluster
487  std::vector<int> jClustCells; //cells with membership in j cluster
488  int overlap = 0; //number of cells appearing in both clusters
489  for (unsigned int k=0; k<fX.size(); k++){
490  if (fU[i][k] >= fParams.Memb()) iClustCells.push_back(k);
491  if (fU[j][k] >= fParams.Memb()) jClustCells.push_back(k);
492  if ((fU[i][k] >= fParams.Memb())&&(fU[j][k] >= fParams.Memb())) ++overlap;
493  }
494 
495  //now computer fraction of smaller cluster that belongs to bigger cluster
496  double frac = 0.0;
497  if (iClustCells.size() <= jClustCells.size()){
498  frac = (double)overlap/iClustCells.size();
499 
500  //if small cluster is above minimum membership in large cluster, merge
501  if (frac >= fParams.MergeMemb()){
502  //add cells to larger cluster with minimum membership
503  for (unsigned int k=0; k<fX.size(); k++){
504  if ((fU[i][k] >= fParams.Memb())&&(fU[j][k] < fParams.Memb())) fU[j][k] = fParams.Memb();
505  }
506  //recompute mean of larger cluster
507  double s = 0.0;
508  double w = 0.0;
509  for (unsigned int k=0; k<fX.size(); ++k) {
510  s += pow(fU[j][k],fParams.Fuzz())*fX[k]/(fdX[k]*fdX[k]);
511  w += pow(fU[j][k],fParams.Fuzz())/(fdX[k]*fdX[k]);
512  }
513  fA[j] = s/w;
514  //now remove the small cluster from the list
515  fA.erase(fA.begin()+i);
516  fU.erase(fU.begin()+i);
517  fD.erase(fD.begin()+i);
518  break;
519  }
520 
521  }
522  // else if the j cluster is the smaller one
523  else{
524  frac = (double)overlap/jClustCells.size();
525 
526  //if small cluster is above minimum membership in large cluster, merge
527  if (frac >= fParams.MergeMemb()){
528  //add cells to larger cluster with minimum membership
529  for (unsigned int k=0; k<fX.size(); k++){
530  if ((fU[j][k] >= fParams.Memb())&&(fU[i][k] < fParams.Memb())) fU[i][k] = fParams.Memb();
531  }
532  //recompute mean of larger cluster
533  double s = 0.0;
534  double w = 0.0;
535  for (unsigned int k=0; k<fX.size(); ++k) {
536  s += pow(fU[i][k],fParams.Fuzz())*fX[k]/(fdX[k]*fdX[k]);
537  w += pow(fU[i][k],fParams.Fuzz())/(fdX[k]*fdX[k]);
538  }
539  fA[i] = s/w;
540  //now remove the small cluster from the list
541  fA.erase(fA.begin()+j);
542  fU.erase(fU.begin()+j);
543  fD.erase(fD.begin()+j);
544  continue;
545  }
546 
547  } //end else
548  }//end if clusters are close in angle
549  }//end loop on j
550 
551  //if a merge was made, continue to see if other merges are possible with that cluster
552  if (size != fA.size()){
553  size = fA.size();
554  continue;
555  }
556 
557  ++i; //if no merge has been made increase i iterator to next cluster
558  }//end while loop on i
559 
560  }
std::vector< std::vector< double > > fD
Dist clust i to sub. j.
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
Atom< double > Fuzz
Definition: FuzzyKMeanAlg.h:32
constexpr T pow(T x)
Definition: pow.h:75
std::vector< double > fdX
Uncertainty in angle.
#define M_PI
Definition: SbMath.h:34
std::vector< double > fX
Angle to vertex.
const XML_Char * s
Definition: expat.h:262
std::vector< double > fA
Cluster centers.
std::vector< std::vector< double > > fU
Assign cluster i to sub. j.
Atom< double > MergeAng
Definition: FuzzyKMeanAlg.h:46
const double j
Definition: BetheBloch.cxx:29
double frac(double x)
Fractional part.
FuzzyKMeanParams fParams
Atom< double > Memb
Definition: FuzzyKMeanAlg.h:35
Float_t w
Definition: plot.C:20
Atom< double > MergeMemb
Definition: FuzzyKMeanAlg.h:47
void fuzz::FuzzyKMeanAlg::SeedA ( unsigned int  n)

A function to calculate the seeds used for the first round of clustering.

Function to produce n seeds for the candidate clusters.

Parameters
nNumber of seeds to calculate

Definition at line 164 of file FuzzyKMeanAlg.cxx.

References angle, fuzz::FuzzyKMeanParams::Bins, release_diff::diff, stan::math::exp(), fA, fdX, fLMax, fParams, fSeed, fX, MECModelEnuComparisons::i, calib::j, M_PI, make_pair(), getGoodRuns4SAM::n, site_stats_from_log::reverse, moon_position_table_new3::second, fillBadChanDBTables::step, and w.

Referenced by DoClustering().

165  {
166  //produce intial seeds evenly distributed about circle
167  //this is done in case there are not enough dense regions
168  //to provide the desired number of seeds
169  double step = M_PI/n;
170  double angle = 0.0;
171  for (unsigned int i=0; i<n; ++i) {
172  fA[i] = angle;
173  angle += step;
174  if (angle > M_PI) angle -= 2*M_PI;
175  }
176 
177  //calculate the angular density in 360 bins from -M_PI to M_PI
178  std::vector<double> w(fParams.Bins());
179  double diff;
180  step = (double)(2.0*M_PI)/fParams.Bins();
181  for (unsigned int k=0; k<w.size(); k++){
182  w[k] = 0.0;
183  for (unsigned int j=0; j<fX.size(); j++){
184  diff = -M_PI + k*step - fX[j];
185  while (diff <-M_PI) diff += 2*M_PI;
186  while (diff > M_PI) diff -= 2*M_PI;
187  w[k] += exp(-diff*diff/(fdX[j]*fdX[j]));
188  }
189  }
190 
191  //now find local maxima in weight distribution
192  fLMax.clear();
193  for (unsigned int k=0; k<w.size(); k++) {
194  if (k==0){
195  if ((w[k] >= w.back()) && (w[k] > w[k+1])) {
196  fLMax.push_back(std::make_pair(w[k],-M_PI+k*step));
197  }
198  }
199  else if (k==w.size()-1) {
200  if ((w[k] >= w[k-1]) && (w[k] > w.front())) {
201  fLMax.push_back(std::make_pair(w[k],-M_PI+k*step));
202  }
203  }
204  else {
205  if ((w[k] >= w[k-1]) && (w[k] > w[k+1])) {
206  fLMax.push_back(std::make_pair(w[k],-M_PI+k*step));
207  }
208  }
209  }
210  // get sorted density list
211  std::sort(fLMax.begin(),fLMax.end());
212  std::reverse(fLMax.begin(),fLMax.end());
213 
214  //populate cluster seeds starting with the maximum density
215  for (unsigned int k=0; k<n; k++){
216  if (k<fLMax.size()) fA[k] = fLMax[k].second;
217  fSeed.push_back(fA[k]);
218  }
219 
220  }
Double_t angle
Definition: plot.C:86
std::vector< double > fSeed
List of used seeds.
std::vector< double > fdX
Uncertainty in angle.
std::pair< Spectrum *, CheatDecomp * > make_pair(SpectrumLoaderBase &loader_data, SpectrumLoaderBase &loader_mc, HistAxis *axis, Cut *cut, const SystShifts &shift, const Var &wei)
Definition: DataMCLoad.C:336
#define M_PI
Definition: SbMath.h:34
std::vector< double > fX
Angle to vertex.
std::vector< std::pair< double, double > > fLMax
List of max densities.
std::vector< double > fA
Cluster centers.
fvar< T > exp(const fvar< T > &x)
Definition: exp.hpp:10
const double j
Definition: BetheBloch.cxx:29
FuzzyKMeanParams fParams
Float_t w
Definition: plot.C:20

Member Data Documentation

std::vector<double> fuzz::FuzzyKMeanAlg::fA
std::vector< std::vector<double> > fuzz::FuzzyKMeanAlg::fD

Dist clust i to sub. j.

Definition at line 131 of file FuzzyKMeanAlg.h.

Referenced by ComputeDij(), ComputeUij(), MakeClusters(), and MergeClusters().

std::vector<double> fuzz::FuzzyKMeanAlg::fDist

Distance from vertex.

Definition at line 119 of file FuzzyKMeanAlg.h.

Referenced by fuzz::FuzzyKVertex::AddProng(), and MakeAngles().

std::vector<double> fuzz::FuzzyKMeanAlg::fdX

Uncertainty in angle.

Definition at line 128 of file FuzzyKMeanAlg.h.

Referenced by ComputeA(), ComputeDij(), FindUnclustered(), MakeAngles(), MergeClusters(), and SeedA().

rb::HitList fuzz::FuzzyKMeanAlg::fHL

List of hits in slice.

Definition at line 116 of file FuzzyKMeanAlg.h.

Referenced by LoadCluster(), and MakeAngles().

std::vector<std::pair<double,double> > fuzz::FuzzyKMeanAlg::fLMax

List of max densities.

Definition at line 146 of file FuzzyKMeanAlg.h.

Referenced by CheckUnique(), FindUnclustered(), and SeedA().

FuzzyKMeanParams fuzz::FuzzyKMeanAlg::fParams
private
std::vector<double> fuzz::FuzzyKMeanAlg::fSeed

List of used seeds.

Definition at line 140 of file FuzzyKMeanAlg.h.

Referenced by CheckUnique(), DoClustering(), FindUnclustered(), and SeedA().

std::vector< std::vector<double> > fuzz::FuzzyKMeanAlg::fU

Assign cluster i to sub. j.

Definition at line 134 of file FuzzyKMeanAlg.h.

Referenced by fuzz::FuzzyKVertex::AddProng(), ComputeA(), ComputeUij(), DoClustering(), FindUnclustered(), Jpfuz(), MakeClusters(), and MergeClusters().

std::vector<std::pair<double,double> > fuzz::FuzzyKMeanAlg::fUnused

Unused cell pos., uncert.

Definition at line 143 of file FuzzyKMeanAlg.h.

Referenced by FindUnclustered().

std::vector<double> fuzz::FuzzyKMeanAlg::fW

Hit weight, currently PE.

Definition at line 122 of file FuzzyKMeanAlg.h.

Referenced by MakeAngles().

std::vector<double> fuzz::FuzzyKMeanAlg::fX

The documentation for this class was generated from the following files: