HitTower.cxx

Go to the documentation of this file.

////////////////////////////////////////////////////////////////////////
/// \file    HitTower.cxx
/// \brief   Plot hits (Raw, Calib, ...) as a tower in 3D
/// \author  messier@indiana.edu
/// \version $Id: HitTower.cxx,v 1.5 2011-11-02 18:56:08 gsdavies Exp $
///////////////////////////////////////////////////////////////////////
#include <algorithm>
#include <cmath>
#include <iostream>

#include "TPolyLine3D.h"

#include "art/Framework/Services/Registry/ServiceHandle.h"

#include "EventDisplay/HitTower.h"
#include "nuevdb/EventDisplayBase/View3D.h"
#include "Geometry/Geometry.h"

namespace evd
{
  //.....................................................................
  bool comp_hitdata(const HitData& h1, const HitData& h2) 
  {
    if (h1.fView  < h2.fView)  return true;
    if (h1.fView  > h2.fView)  return false;
    if (h1.fPlane < h2.fPlane) return true;
    if (h1.fPlane > h2.fPlane) return false;
    if (h1.fCell  < h2.fCell)  return true;
    if (h1.fCell  > h2.fCell)  return false;
    return h1.fT < h2.fT;
  }
  
  //......................................................................
  HitData::HitData(int v, int p, int c, float tp, float zp, float q, float t) :
    fView(v),
    fPlane(p),
    fCell(c),
    fTpos(tp),
    fZpos(zp),
    fQ(q),
    fT(t)
  {}

  //......................................................................

  void HitTower::AddHit(int v, int p, int c, float tp, float zp, 
                        float q, float t) 
  {
    fHitData.push_back(HitData(v,p,c,tp,zp,q,t));
  }

  //......................................................................

  void HitTower::Draw(evdb::View3D* view)
  {
    unsigned int i;

    art::ServiceHandle<geo::Geometry> geom;

    // Sort hits by view, plane, cell, time
    std::sort(fHitData.begin(), fHitData.end(), comp_hitdata);
  
    //                                  _
    // Render hits in this style:     _| |_
    //                              _|     |_
    //
    TPolyLine3D* polyline = 0;
    int          ipoint   = 0;
    for (i=0; i<fHitData.size(); ++i) {
      // Get information about the current cell in the list
      const HitData& hit0    = fHitData[i];
      double         pos0[3] = {0,0,0};
      double         pos1[3] = {0,0,0};
      double         dpos[3] = {0,0,0};
      geo::View_t    v       = geo::kXorY;
      geom->CellInfo(hit0.fPlane, hit0.fCell, &v, pos0, dpos);
      if (hit0.fView == geo::kX) {
        pos0[0] = hit0.fTpos - dpos[0];
        pos0[1] = geom->DetHalfHeight()+fQscale*hit0.fQ;
        pos0[2] = hit0.fZpos;

        pos1[0] = hit0.fTpos + dpos[0];
        pos1[1] = pos0[1];
        pos1[2] = pos0[2];
      }
      else if (hit0.fView == geo::kY) {
        pos0[0] = geom->DetHalfWidth()+fQscale*hit0.fQ;
        pos0[1] = hit0.fTpos-dpos[1];
        pos0[2] = hit0.fZpos;

        pos1[0] = pos0[0];
        pos1[1] = hit0.fTpos+dpos[1];
        pos1[2] = pos0[2];
      }
      else abort();
    
      // Get information about the previous hit in the list
      bool connectL = false;
      if (i>0) {
        const HitData& hitL = fHitData[i-1];
        double tcut = 9E9;
        bool samePlane = (hit0.fPlane==hitL.fPlane);
        bool adjaCell  = (abs(hit0.fCell-hitL.fCell)==1);
        bool sameTime  = (fabs(hit0.fT-hitL.fT)<tcut);
        if (samePlane && adjaCell && sameTime) connectL = true;
      }

      // Get information about the next hit in the list
      bool connectR = false;
      if (i+1<fHitData.size()) {
        const HitData& hitR = fHitData[i+1];
        double tcut = 9E9;
        bool samePlane = (hit0.fPlane==hitR.fPlane);
        bool adjaCell  = (abs(hit0.fCell-hitR.fCell)==1);
        bool sameTime  = (fabs(hit0.fT-hitR.fT)<tcut);
        if (samePlane && adjaCell && sameTime) connectR = true;
      }

      // If we don't have a left connection, then we need to start a new
      // line
      if (connectL == false) {
        // How long a lead in stub to draw? 2.0 would fill the whole
        // adjacent cell
        double stubfrac = 1.90;

        // The coordinates for the two points on the stub
        double x0, y0, z0;
        double x1, y1, z1;
        if (hit0.fView == geo::kX) {
          x1 = pos0[0];
          y1 = geom->DetHalfHeight();
          z1 = pos0[2];
        
          x0 = x1 - stubfrac*dpos[0];
          y0 = y1;
          z0 = z1;
        }
        else if (hit0.fView == geo::kY) {
          x1 = geom->DetHalfWidth();
          y1 = pos0[1];
          z1 = pos0[2];
        
          x0 = x1;
          y0 = y1 - stubfrac*dpos[1];
          z0 = z1;
        }
        else abort();
      
        // Make a new line and add two points to the left of the first hit
        ipoint = 0;
        int c = kRed-9;
        int w = 1;
        int s = 1;
        polyline = &view->AddPolyLine3D(0,c,w,s);
        polyline->SetPoint(ipoint,x0,y0,z0); ++ipoint;
        polyline->SetPoint(ipoint,x1,y1,z1); ++ipoint;
      }

      // At this point in the code we better have a polyline or
      // something has gone wrong...
      if (polyline==0) {
        std::cerr << __FILE__ ":" << __LINE__ 
                  << " Polyline=0? connectL=" << connectL << std::endl;
      }

      // Always add two points to the line for the current hit
      polyline->SetPoint(ipoint, pos0[0], pos0[1], pos0[2]); ++ipoint;
      polyline->SetPoint(ipoint, pos1[0], pos1[1], pos1[2]); ++ipoint;
    
      // If we don't have a right connection, then we need to end this
      // line
      if (connectR == false) {
        // How long a lead in stub to draw? 2 would fill the whole
        // adjact cell
        double stubfrac = 1.90;
      
        // The coordinates on the two ends of the stub
        double x0, y0, z0;
        double x1, y1, z1;
        if (hit0.fView == geo::kX) {
          x0 = pos1[0];
          y0 = geom->DetHalfHeight();
          z0 = pos1[2];
        
          x1 = x0 + stubfrac*dpos[0];
          y1 = y0;
          z1 = z0;
        }
        else if (hit0.fView == geo::kY) {
          x0 = geom->DetHalfWidth();
          y0 = pos1[1];
          z0 = pos1[2];
        
          x1 = x0;
          y1 = y0 + stubfrac*dpos[1];
          z1 = z0;
        }
        else abort();
        polyline->SetPoint(ipoint, x0, y0, z0); ++ipoint;
        polyline->SetPoint(ipoint, x1, y1, z1); ++ipoint;
      }
    } // loop on i in fHitData
  }
} // end namespace evd
////////////////////////////////////////////////////////////////////////