SimpleClustering.cxx

Go to the documentation of this file.

#include "SimpleClustering.h"

tss::Cluster2D::Cluster2D(const std::vector<const tss::Hit2D*>& hits)
  : fDenseStart(false), fDenseEnd(false), fIsEM(false)
{
  fHits.reserve(hits.size());
  for (size_t h = 0; h < hits.size(); ++h)
    fHits.push_back(hits[h]);
}
// ------------------------------------------------------

const tss::Hit2D* tss::Cluster2D::release_at(size_t idx)
{
  const tss::Hit2D* hit = 0;
  if (idx < fHits.size()) {
    hit = fHits[idx];
    fHits.erase(fHits.begin() + idx);
  }
  return hit;
}
// ------------------------------------------------------

bool tss::Cluster2D::release(const tss::Hit2D* hit)
{
  for (size_t h = 0; h < fHits.size(); ++h)
    if (fHits[h] == hit) {
      fHits.erase(fHits.begin() + h);
      return true;
    }
  return false;
}

const tss::Hit2D* tss::Cluster2D::closest(const TVector2& p2d, size_t& idx) const
{
  idx = 0;
  if (!fHits.size()) return 0;

  const tss::Hit2D* hout = fHits.front();
  double d, dmin = pma::Dist2(hout->Point2D(), p2d);
  for (size_t h = 1; h < fHits.size(); ++h) {
    d = pma::Dist2(fHits[h]->Point2D(), p2d);
    if (d < dmin) {
      dmin = d;
      hout = fHits[h];
      idx = h;
    }
  }
  return hout;
}
// ------------------------------------------------------

const tss::Hit2D* tss::Cluster2D::outermost(size_t& idx) const
{
  idx = 0;
  if (!fHits.size()) return 0;

  TVector2 mean(0., 0.);
  for (size_t h = 0; h < fHits.size(); ++h) {
    mean += fHits[h]->Point2D();
  }
  mean *= 1.0 / fHits.size();

  const tss::Hit2D* hout = fHits.front();
  double d, dmax = pma::Dist2(hout->Point2D(), mean);
  for (size_t h = 1; h < fHits.size(); ++h) {
    d = pma::Dist2(fHits[h]->Point2D(), mean);
    if (d > dmax) {
      dmax = d;
      hout = fHits[h];
      idx = h;
    }
  }
  return hout;
}
// ------------------------------------------------------

const TVector2 tss::Cluster2D::min(void) const
{

  TVector2 minimum = fHits[0]->Point2D();

  for (size_t i = 1; i < size(); ++i) {
    const TVector2 h = fHits[i]->Point2D();

    if (h.X() < minimum.X()) minimum.Set(h.X(), h.Y());
    if (h.Y() < minimum.Y()) minimum.Set(minimum.X(), h.Y());
  }

  return minimum;
}

// ------------------------------------------------------

const TVector2 tss::Cluster2D::max(void) const
{

  TVector2 maximum = fHits[0]->Point2D();

  for (size_t i = 1; i < size(); ++i) {
    const TVector2 h = fHits[i]->Point2D();

    if (h.X() > maximum.X()) maximum.Set(h.X(), h.Y());
    if (h.Y() > maximum.Y()) maximum.Set(maximum.X(), h.Y());
  }

  return maximum;
}

// ------------------------------------------------------

bool tss::Cluster2D::has(const tss::Hit2D* hit) const
{
  for (size_t i = 0; i < fHits.size(); ++i)
    if (fHits[i] == hit) return true;
  return false;
}
// ------------------------------------------------------

double tss::Cluster2D::dist2(const TVector2& p2d) const
{
  if (fHits.size()) {
    double d2, min_d2 = pma::Dist2(fHits.front()->Point2D(), p2d);
    for (size_t i = 1; i < fHits.size(); ++i) {
      d2 = pma::Dist2(fHits[i]->Point2D(), p2d);
      if (d2 < min_d2) { min_d2 = d2; }
    }
    return min_d2;
  }
  else
    return 0.0;
}
// ------------------------------------------------------

double tss::Cluster2D::dist2(const TVector2& p2d, size_t& hIdx) const
{
  hIdx = 0;
  if (fHits.size()) {
    double d2, min_d2 = pma::Dist2(fHits.front()->Point2D(), p2d);
    for (size_t i = 1; i < fHits.size(); ++i) {
      d2 = pma::Dist2(fHits[i]->Point2D(), p2d);
      if (d2 < min_d2) {
        min_d2 = d2;
        hIdx = i;
      }
    }
    return min_d2;
  }
  else
    return 0.0;
}
// ------------------------------------------------------

double tss::Cluster2D::dist2(const tss::Cluster2D& clu) const
{
  if (fHits.size()) {
    double d2, min_d2 = clu.dist2(fHits.front()->Point2D());
    for (size_t i = 1; i < fHits.size(); ++i) {
      d2 = clu.dist2(fHits[i]->Point2D());
      if (d2 < min_d2) { min_d2 = d2; }
    }
    return min_d2;
  }
  else
    return 0.0;
}
// ------------------------------------------------------

// ------------------------------------------------------
// ------------------------------------------------------
// ------------------------------------------------------

bool tss::SimpleClustering::hitsTouching(const tss::Hit2D& h1, const tss::Hit2D& h2) const
{
  if ((h1.Wire() == h2.Wire()) && (h1.PeakTime() == h2.PeakTime())) return false;

  bool touches = false;
  if ((h1.Wire() >= h2.Wire() - 1) && (h1.Wire() <= h2.Wire() + 1)) {
    if (((h2.StartTick() <= h1.StartTick()) && (h1.StartTick() <= h2.EndTick() + 1)) ||
        ((h2.StartTick() <= h1.EndTick() + 1) && (h1.EndTick() <= h2.EndTick())) ||
        ((h2.StartTick() >= h1.StartTick()) && (h1.EndTick() >= h2.EndTick()))) {
      touches = true;
    }
  }
  return touches;
}
// ------------------------------------------------------

bool tss::SimpleClustering::hitsTouching(const tss::Cluster2D& c1, const tss::Hit2D& h2) const
{
  for (size_t i = 0; i < c1.size(); i++) {
    if (hitsTouching(c1[i], h2)) return true;
  }
  return false;
}
// ------------------------------------------------------

bool tss::SimpleClustering::hitsTouching(const tss::Cluster2D& c1, const tss::Cluster2D& c2) const
{
  for (unsigned int i = 0; i < c1.size(); i++) {
    if (hitsTouching(c2, c1[i])) return true;
  }
  return false;
}
// ------------------------------------------------------

void tss::SimpleClustering::merge(std::vector<tss::Cluster2D>& clusters) const
{
  bool merged = true;
  while (merged) {
    merged = false;

    size_t i = 0;
    while (i < clusters.size() - 1) {
      size_t j = i + 1;
      while (j < clusters.size()) {
        if (hitsTouching(clusters[i], clusters[j])) {
          clusters[i].hits().reserve(clusters[i].size() + clusters[i].size());
          for (size_t h = 0; h < clusters[j].size(); ++h)
            clusters[i].hits().push_back(clusters[j].hits()[h]);
          clusters.erase(clusters.begin() + j);
          merged = true;
        }
        else
          ++j;
      }
      ++i;
    }
  }
}
// ------------------------------------------------------

std::vector<tss::Cluster2D> tss::SimpleClustering::run(const std::vector<tss::Hit2D>& inp) const
{
  std::vector<tss::Cluster2D> result;
  for (size_t h = 0; h < inp.size(); ++h) {
    bool found = false;
    for (size_t r = 0; r < result.size(); ++r)
      if (hitsTouching(result[r], inp[h])) {
        result[r].hits().push_back(&(inp[h]));
        found = true;
        break;
      }
    if (!found) {
      result.push_back(tss::Cluster2D());
      result.back().hits().push_back(&(inp[h]));
    }
  }
  merge(result);

  return result;
}
// ------------------------------------------------------

std::vector<tss::Cluster2D> tss::SimpleClustering::run(const tss::Cluster2D& inp) const
{
  std::vector<tss::Cluster2D> result;
  for (size_t h = 0; h < inp.size(); ++h) {
    bool found = false;
    for (size_t r = 0; r < result.size(); ++r)
      if (hitsTouching(result[r], inp[h])) {
        result[r].hits().push_back(inp.hits()[h]);
        found = true;
        break;
      }
    if (!found) {
      result.push_back(tss::Cluster2D());
      result.back().hits().push_back(inp.hits()[h]);
    }
  }
  merge(result);

  return result;
}
// ------------------------------------------------------