RawDataDrawer.cxx

Go to the documentation of this file.

///
/// \file    RawDataDrawer.cxx
/// \brief   Class to aid in the rendering of RawData objects
/// \author  messier@indiana.edu
/// \version $Id: RawDataDrawer.cxx,v 1.30 2012-08-28 22:00:57 denis Exp $
///
#include <cmath>

#include "TBox.h"
#include "TLine.h"
#include "TText.h"
#include "TH1F.h"
#include "TMarker3DBox.h"
#include "TPolyLine3D.h"
#include "TPolyMarker3D.h"

#include "art/Framework/Principal/Event.h"
#include "art/Framework/Services/Registry/ServiceHandle.h"
#include "canvas/Persistency/Common/Ptr.h"

#include "nuevdb/EventDisplayBase/Colors.h"
#include "nuevdb/EventDisplayBase/ColorScale.h"

#include "DAQChannelMap/HardwareDisplay.h"
#include "ChannelInfo/BadChanList.h"
#include "CMap/service/CMapService.h"
#include "EventDisplay/HitTower.h"
#include "EventDisplay/RawDataDrawer.h"
#include "EventDisplayServices/RawDrawingOptions.h"
#include "EventDisplayServices/GeometryDrawingOptions.h"
#include "EventDisplay/GeoTransform.h"
#include "nuevdb/EventDisplayBase/View2D.h"
#include "nuevdb/EventDisplayBase/View3D.h"
#include "Geometry/Geometry.h"
#include "GeometryObjects/CellGeo.h"
#include "RawData/RawDigit.h"
#include "RecoBase/CellHit.h"
#include "RawData/RawTrigger.h"
#include "Calibrator/Calibrator.h"
#include "EventDisplayServices/SliceNavigator.h"

namespace {
  // Utility function to make uniform error messages.
  void writeErrMsg(const char* fcn,
                   cet::exception const& e)
  {
    std::cout << "RawDataDrawer::" << fcn
              << " failed with message:\n"
              << e
              << std::endl;
  }
}

namespace evd
{
  //
  // Conversion constant from NOvA clock ticks to micro-seconds
  //
  static double kUSEC_PER_TDC = 1.e6/64000000.0;

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

  RawDataDrawer::RawDataDrawer() { }

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

  RawDataDrawer::~RawDataDrawer() { }

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

  void RawDataDrawer::FillTQHisto(const art::Event& evt,
                                  TH1F*             thisto,
                                  TH1F*             qhisto) 
  {
    // Check if we're supposed to draw raw hits at all
    static art::ServiceHandle<evd::RawDrawingOptions> drawopt;
    if (drawopt->fWhichHits!=evd::RawDrawingOptions::kRAW_HITS) {
      return;
    }

    std::vector<const rawdata::RawDigit*> rawhits;
    this->GetRawDigitsAllLabels(evt, rawhits);
      
    for (unsigned int i=0; i<rawhits.size(); ++i) {
      double t = rawhits[i]->TDC() * kUSEC_PER_TDC;
      double q = rawhits[i]->ADC();
      if (thisto!=0) thisto->Fill(t);
      if (qhisto!=0) {
        //
        // Only histogram ADC counts that are in time with the
        // selected time range
        //
        if (t>drawopt->fTimeRange[0] && t<drawopt->fTimeRange[1]) {
          qhisto->Fill(q);
        }
      }
    }
  }

  //......................................................................
  ///
  /// Draw raw data according to hardware address
  ///
  void RawDataDrawer::RawDigitHardware(const art::Event& evt,
                                       evdb::View2D*     view)
  {
    art::ServiceHandle<evd::RawDrawingOptions> drawopt;
    if (drawopt->fWhichHits!=evd::RawDrawingOptions::kRAW_HITS) return;

    art::ServiceHandle<geo::Geometry>               geo;
    art::ServiceHandle<evd::GeometryDrawingOptions> geoopt;
    art::ServiceHandle<evd::SliceNavigator>         nav;

    daqchannelmap::HardwareDisplay hd;
    hd.SetupDet(geo->DetId());
    
    std::vector<const rb::Cluster*> slicesToNavigate = nav->SlicesToNavigate(evt);

    // channel, tdc
    std::set<std::pair<uint32_t, int16_t>> unfiltHits;

    int n = slicesToNavigate.size();
    for(int i=0 ; i < n; ++i ){
      // If DimOthers option is chosen, hits from slices
      // other than the current slice should be dimmed.

      if( nav->DimOthers() && !nav->FullEventMode() ){
        if( i != nav->CurrentSliceIndex() ) 
          continue;
      }
      for( int c = 0, nc = slicesToNavigate[i]->NCell(); c < nc; ++c){
        rb::CellHit chit = *(slicesToNavigate[i]->Cell(c));

        unfiltHits.emplace(chit.Channel(), chit.TDC());
      }
    }

    unsigned int ix; // hardware x coordinate
    unsigned int iy; // hardware y coordinate

    std::vector< const rawdata::RawDigit*> rawhits;
    this->GetRawDigitsAllLabels(evt, rawhits);
      
    for (unsigned int i=0; i<rawhits.size(); ++i) {
      unsigned int ch = rawhits[i]->DaqChannel();
      unsigned int diblock = (ch&0x0FC00000)>>22;
      unsigned int dcm     = (ch&0x003F0000)>>16;
      unsigned int feb     = (ch&0x0000FF00)>>8;
      unsigned int pix     = (ch&0x000000FF)>>0;
        
      hd.DAQToXY(diblock, dcm, feb, pix, &ix, &iy);
      double t = 0;
      double q = 0;
        
      t = rawhits[i]->TDC() * kUSEC_PER_TDC;
      q = rawhits[i]->ADC();
      if (q<=0.0) q=1.;
        
      //
      // Draw markers representing hits
      //
      int ghosted;
      int co = this->HitColor(rawhits[i], unfiltHits, &ghosted);
      if (drawopt->fRawDrawingOpt &
          evd::RawDrawingOptions::kSUPPRESS_GHOSTED_HITS) {
        if (ghosted) continue;
      }
      double sf = this->ScaleFactor(t, q);

      TBox& b1 = view->AddBox(ix-0.5*sf,iy-0.5*sf,ix+0.5*sf,iy+0.5*sf);
      b1.SetFillStyle(1001);
      b1.SetFillColor(co);    
      b1.SetBit(kCannotPick);

      if (geoopt->fOutline&evd::GeometryDrawingOptions::kOUTLINE_CELLS) {
        TBox& b2 = view->AddBox(ix-0.5*sf,iy-0.5*sf,ix+0.5*sf,iy+0.5*sf);
        b2.SetFillStyle(0);
        if (ghosted) b2.SetLineColor(17); // Greyed out hit
        else         b2.SetLineColor(14);
        b2.SetBit(kCannotPick);
      }
    }// end loop over hits
    
    // Now draw hardware boundaries and labels
    
    // Vertical lines marking the FEB boundaries. Also mark diblock
    // boundaries periodically
    
    for (unsigned int i=0; i<=hd.fNdiblocks*hd.fNfebCols; ++i) {
      int lw = 1;
      int co = 33;
      if (i%hd.fNfebCols==0) { lw = 1; co = 38; }
      TLine& l = view->AddLine(i*hd.fNpixCols-0.5, -0.5,
                               i*hd.fNpixCols-0.5, hd.fYmax-0.5);
      l.SetLineStyle(1);
      l.SetLineColor(co);
      l.SetLineWidth(lw);
      l.SetBit(kCannotPick);
    }
    for (unsigned int i=0; i<(hd.fNdcmTop+hd.fNdcmSide)*hd.fNfebRows; ++i) {
      int lw = 1;
      int co = 33;
      if (i%hd.fNfebRows==0) { lw = 1; co = 38; }
      TLine& l = view->AddLine(-0.5,        i*hd.fNpixRows-0.5,
                               hd.fXmax-0.5,i*hd.fNpixRows-0.5);
      l.SetLineStyle(1);
      l.SetLineColor(co);
      l.SetLineWidth(lw);
      l.SetBit(kCannotPick);
    }

    for (unsigned int i=1; i<=hd.fNdiblocks; ++i) {
      char buff[128];
      sprintf(buff,"DB%.2u",i);
      iy = -hd.fYmax/200;
      hd.DAQToXY(i,1,0,0,&ix,&iy);
      TText& t = view->AddText(ix,-1,buff);
      t.SetTextSize(0.04);
      t.SetTextAlign(13);
    }

    for (unsigned int i=1; i<=hd.fNdcmTop+hd.fNdcmSide; ++i) {
      char buff[128];
      sprintf(buff,"DCM%.2u",i);
      hd.DAQToXY(1,i,32,0,&ix,&iy);
      ix = -hd.fXmax/200;
      TText& t = view->AddText(-1,iy,buff);
      t.SetTextSize(0.04);
      t.SetTextAlign(31);
    }
  }

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

  void RawDataDrawer::RawDigit2D(const art::Event& evt,
                                 evdb::View2D*     xzview,
                                 evdb::View2D*     yzview) 
  {
    // Check if we're supposed to draw raw hits at all
    art::ServiceHandle<evd::RawDrawingOptions> drawopt;
    if (drawopt->fWhichHits!=evd::RawDrawingOptions::kRAW_HITS) return;

    art::ServiceHandle<geo::Geometry>               geo;
    art::ServiceHandle<cmap::CMap>                  cmap;
    art::ServiceHandle<evd::GeometryDrawingOptions> geoopt;
    art::ServiceHandle<evd::SliceNavigator> nav;

    std::vector<const rb::Cluster*> slicesToNavigate = nav->SlicesToNavigate(evt);
    std::set<std::pair<uint32_t, int16_t>> unfiltHits;

    for(int i=0, n = slicesToNavigate.size(); i < n; ++i ){
      // If DimOthers option is chosen, hits from slices
      // other than the current slice should be dimmed.
      if( nav->DimOthers() && !nav->FullEventMode() ){
        if( i != nav->CurrentSliceIndex() ) 
          continue;
      }
      for( int c = 0, nc = slicesToNavigate[i]->NCell(); c < nc; ++c){
        rb::CellHit chit = *(slicesToNavigate[i]->Cell(c));
        unfiltHits.emplace(chit.Channel(), chit.TDC());
      }
    }

    fXmin =  9.E9;
    fYmin =  9.E9;
    fZmin =  9.E9;
    fXmax = -9.E9;
    fYmax = -9.E9;
    fZmax = -9.E9;

    // Draw the dimmed hits first to ensure selected hits are above them
    for(bool drawDimmed: {true, false}){
      std::vector< const rawdata::RawDigit*> rawhits;
      this->GetRawDigitsAllLabels(evt, rawhits);
      
      for (unsigned int i=0; i<rawhits.size(); ++i) {
        unsigned int plane = cmap->GetPlane(rawhits[i]);
        unsigned int cell  = cmap->GetCell(rawhits[i]);
        
        if (plane>=geo->NPlanes()) {
          std::cerr << "*- Illegal plane number " << plane
                    << " DaqChannel = " 
                    << std::hex << rawhits[i]->DaqChannel() 
                    << std::endl;
          std::cerr << std::dec;
          continue;
        }
        const geo::PlaneGeo *p = geo->Plane(plane);
        
        if (cell>=p->Ncells()) {
          std::cerr << "*- Illegal cell number " << cell 
                    << " in plane " << plane
                    << "*- DaqChannel = " 
                    << std::hex << rawhits[i]->DaqChannel()
                    << std::endl;
          std::cerr << std::dec;
          continue;
        }
        const geo::CellGeo *c = p->Cell(cell);
        
        double xyz[3];
        double tpos = 0.0, zpos = 0.0;
        double dt, dz;
        dz = c->HalfD();
        dt = c->HalfW();
        evdb::View2D* v = 0;
        if (p->View()==geo::kX) {
          c->GetCenter(xyz);
          GeoTransform::XYZ(xyz);
          
          tpos = xyz[0];
          zpos = xyz[2];
      
          if (tpos<fXmin) fXmin = tpos;
          if (tpos>fXmax) fXmax = tpos;
          if (zpos<fZmin) fZmin = zpos;
          if (zpos>fZmax) fZmax = zpos;
        }
        if (p->View()==geo::kY) {
          c->GetCenter(xyz);
          GeoTransform::XYZ(xyz);

          tpos = xyz[1];
          zpos = xyz[2];

          if (tpos<fYmin) fYmin = tpos;
          if (tpos>fYmax) fYmax = tpos;
          if (zpos<fZmin) fZmin = zpos;
          if (zpos>fZmax) fZmax = zpos;
        }
        // Decide if this hit should be drawn on this view
        if (p->View()==geo::kX && xzview==0) continue;
        if (p->View()==geo::kY && yzview==0) continue;
        if (xzview!=0) v = xzview;
        if (yzview!=0) v = yzview;
    
        double t = 0;
        double q = 0;
        t = rawhits[i]->TDC()*kUSEC_PER_TDC;
        q = rawhits[i]->ADC();

        // Hack for now...
        if (q<=0.0) q=100.;
        
        //
        // Draw boxes representing hit cells
        //
        if(DimFiltered(rawhits[i], unfiltHits) != drawDimmed) continue;

        int ghosted;
        int co = this->HitColor( rawhits[i], unfiltHits, &ghosted);
        if (drawopt->fRawDrawingOpt & 
            evd::RawDrawingOptions::kSUPPRESS_GHOSTED_HITS) {
          if (ghosted) continue;
        }

        double sf = this->ScaleFactor(t, q);
        TBox& b1 = v->AddBox(zpos-sf*dz,tpos-sf*dt,zpos+sf*dz,tpos+sf*dt);
        b1.SetFillStyle(1001);
        b1.SetFillColor(co);    
        b1.SetBit(kCannotPick);

        if (geoopt->fOutline&evd::GeometryDrawingOptions::kOUTLINE_CELLS) {
          TBox& b2 = v->AddBox(zpos-dz,tpos-dt,zpos+dz,tpos+dt);
          b2.SetFillStyle(0);
          if (ghosted) b2.SetLineColor(17); // Greyed out hit
          else         b2.SetLineColor(14);
          b2.SetBit(kCannotPick);
        }
      }//end loop over raw hits
    } // end for drawDimmed
  }

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

  void RawDataDrawer::
  Draw3DHitBoxes(std::vector<const rawdata::RawDigit*>& rawhits,
                 const std::set<std::pair<uint32_t, int16_t>>& unfiltHits,
                 evdb::View3D*                          view)
  {
    art::ServiceHandle<geo::Geometry> geom;
    art::ServiceHandle<cmap::CMap>    cmap;
    for (unsigned int i=0; i<rawhits.size(); ++i) {
      double t = 0;
      double q = 0;
      t = rawhits[i]->TDC()*kUSEC_PER_TDC;
      q = rawhits[i]->ADC();
      
      // Get the cell geometry for the hit
      int         iplane = cmap->GetPlane(rawhits[i]);
      int         icell  = cmap->GetCell(rawhits[i]);
      double      xyz[3];
      double      dxyz[3];
      geo::View_t v;
      geom->CellInfo(iplane, icell, &v, xyz, dxyz);
      
      int ghosted;
      int co = this->HitColor(rawhits[i], unfiltHits, &ghosted);
      if (ghosted) continue;
      
      double sf = this->ScaleFactor(t, q);
      double h  = 0.0;

      double x0, y0, z0;
      double dx, dy, dz;
      if (v==geo::kX) {
        x0 = xyz[0];
        y0 = geom->DetHalfHeight()+0.5*h;
        z0 = xyz[2];
        dx = sf*dxyz[0];
        dy = 0.5*h;
        dz = sf*dxyz[2];
      }      
      else {
        x0 = geom->DetHalfWidth()+0.5*h;
        y0 = xyz[1];
        z0 = xyz[2];
        dx = 0.5*h;
        dy = sf*dxyz[1];
        dz = sf*dxyz[2];
      }      
      TMarker3DBox& m = view->AddMarker3DBox(x0, y0, z0, dx, dy, dz);
      m.SetLineStyle(1);
      m.SetLineColor(co);
    }
  }

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

  void RawDataDrawer::RawDigit3D(const art::Event& evt,
                                 evdb::View3D*     view)
  {
    // Check if we're supposed to draw raw hits at all
    art::ServiceHandle<evd::RawDrawingOptions> drawopt;
    if (drawopt->fWhichHits!=evd::RawDrawingOptions::kRAW_HITS) return;

    art::ServiceHandle<geo::Geometry> geom;
    art::ServiceHandle<evd::SliceNavigator> nav;

    std::vector<const rb::Cluster*> slicesToNavigate = nav->SlicesToNavigate(evt);
    std::set<std::pair<uint32_t, int16_t>> unfiltHits;

    for (int i=0, n = slicesToNavigate.size(); i<n; ++i){
      //
      // If DimOthers option is chosen, hits from slices other than
      // the current slice should be dimmed.
      //
      if (nav->DimOthers() && !nav->FullEventMode()) {
        if (i!=nav->CurrentSliceIndex()) continue;
      }
      for(int c=0, nc=slicesToNavigate[i]->NCell(); c<nc; ++c){
        rb::CellHit chit = *(slicesToNavigate[i]->Cell(c));
        unfiltHits.emplace(chit.Channel(), chit.TDC());
      }
    } // end loop over slices to navigate

    std::vector< const rawdata::RawDigit*> rawhits;
    this->GetRawDigitsAllLabels(evt, rawhits);

    this->Draw3DHitBoxes(rawhits, unfiltHits, view);
    
    //
    // A stick figure for scale
    //
    this->StickFigure(view,
                       geom->DetHalfWidth(),
                      -geom->DetHalfHeight(),
                      0.0);
  }

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

  void RawDataDrawer::GetLimits(double* xmin, double* xmax,
                                double* ymin, double* ymax,
                                double* zmin, double* zmax) 
  {
    *xmin = fXmin;
    *xmax = fXmax;
    *ymin = fYmin;
    *ymax = fYmax;
    *zmin = fZmin;
    *zmax = fZmax;
  }

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

  template<class T> bool RawDataDrawer::
  GetRawDigits(const art::Event& evt,
               const std::string& which,
               const std::string& instance,
               std::vector<const rawdata::RawDigit*>& rawhits)
  {
    rawhits.clear();

    // Temporary vector of rawdigits
    std::vector<const rawdata::RawDigit*> temp;

    art::Handle<std::vector<T>> rdcol;
    evt.getByLabel(which, instance, rdcol);

    if(rdcol.failedToGet()) return false;

    // Extract rawdigits from Handle
    temp.reserve(rdcol->size());
    for(const T& dig: *rdcol) temp.push_back(&dig);
      
    // Fill in the temporary vector of rawdigits to output reference rawhits
    temp.swap(rawhits);

    // check for bad channels
    art::ServiceHandle<evd::RawDrawingOptions> drawopt;
    if (drawopt->fRawDrawingOpt & evd::RawDrawingOptions::kMASK_BAD_CHANNELS) {
      art::ServiceHandle<chaninfo::BadChanList> bcl;
      bcl->Apply(rawhits);
    }

    return true;
  }

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

  bool RawDataDrawer::
  GetRawDigitsMulti(const art::Event& evt,
                    std::vector<const rawdata::RawDigit*>& rawhits)
  {
    rawhits.clear();

    // Temporary vector of rawdigits
    std::vector<const rawdata::RawDigit*> temp;

    // Try to find digits in other folders
    std::vector<art::Handle<std::vector<rawdata::RawDigit>>> rdcol_vec;
    evt.getManyByType(rdcol_vec);

    if(rdcol_vec.empty()) return false;

    for(const art::Handle<std::vector<rawdata::RawDigit>>& digs: rdcol_vec)
      for(const rawdata::RawDigit& dig: *digs)
        temp.push_back(&dig);
      
    /// Fill in the temporary vector of rawdigits to output reference rawhits
    temp.swap(rawhits);

    // check for bad channels
    art::ServiceHandle<evd::RawDrawingOptions> drawopt;
    if (drawopt->fRawDrawingOpt & evd::RawDrawingOptions::kMASK_BAD_CHANNELS) {
      art::ServiceHandle<chaninfo::BadChanList> bcl;
      bcl->Apply(rawhits);
    }

    return true;
  }

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

  void RawDataDrawer::GetRawDigitsAllLabels(const art::Event& evt,
                                            std::vector<const rawdata::RawDigit*>& rawhits)
  {
    static art::ServiceHandle<evd::RawDrawingOptions> drawopt;

    rawhits.clear();

    bool any = false;
    for(unsigned imod = 0; imod < drawopt->fRawDigitsModules.size(); ++imod){
      std::vector<const rawdata::RawDigit*> digs;
      if(GetRawDigits<rawdata::RawDigit>(evt,
                                         drawopt->fRawDigitsModules[imod],
                                         drawopt->fRawDigitsInstances[imod],
                                         digs)){
        any = true;
        rawhits.insert(rawhits.end(), digs.begin(), digs.end());
      }
    }
    if(any) return;

    if(GetRawDigitsMulti(evt, rawhits)) return;

    for(unsigned imod = 0; imod < drawopt->fCellHitsModules.size(); ++imod){
      std::vector<const rawdata::RawDigit*> digs;
      if(GetRawDigits<rb::CellHit>(evt,
                                   drawopt->fCellHitsModules[imod],
                                   drawopt->fCellHitsInstances[imod],
                                   digs)){
        any = true;
        rawhits.insert(rawhits.end(), digs.begin(), digs.end());
      }
    }

    if(!any){
      std::cout << __PRETTY_FUNCTION__
                << ": Could not find any hits. Looked under labels";
      for(const std::string& m: drawopt->fRawDigitsModules) std::cout << " '" << m << "'";
      std::cout << " (for RawDigits) and";
      for(const std::string& m: drawopt->fCellHitsModules) std::cout << " '" << m << "'";
      std::cout << " (for CellHits)" << std::endl;
    }
  }

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

  int RawDataDrawer::GetRawTrigger(const art::Event&                    e,
                                   art::PtrVector<rawdata::RawTrigger>& r)
  {
    r.clear();
    art::PtrVector<rawdata::RawTrigger> temp;

    // using getManyByType so that we don't have to specify a module 
    // label
    std::vector< art::Handle< std::vector<rawdata::RawTrigger> > > rtcol;
  
    try{
      e.getManyByType(rtcol);
      for(size_t rtc = 0; rtc < rtcol.size(); ++rtc){
        for(unsigned int i=0; i < rtcol[rtc]->size(); ++i){
          art::Ptr<rawdata::RawTrigger> rt(rtcol[rtc], i);
          temp.push_back(rt);
        }
      }
      temp.swap(r);
    }// end try
    catch(cet::exception& e) {
      writeErrMsg("GetRawTrigger", e);
    }

    return r.size();
  }

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

  ///
  /// Determine the color to use for the hit. If the hit is flagged
  /// to be ghosted out by either its time or charge, ghost the hit
  ///
  int RawDataDrawer::HitColor(const rawdata::RawDigit* hit,
                              const std::set<std::pair<uint32_t, int16_t>>& unfiltered,
                              int* ghosted)
  {
    static art::ServiceHandle<evd::RawDrawingOptions>  drawopt;
    static art::ServiceHandle<evd::SliceNavigator>     sliceopt;

    int co = 0;
    double t = hit->TDC() * kUSEC_PER_TDC;
    double q = hit->ADC();

    static art::ServiceHandle<evdb::Colors> colors;
    int t_color = colors->Scale("RawT").GetColor(t);
    int q_color = colors->Scale("RawQ").GetColor(q);

    *ghosted =
      !(colors->Scale("RawT").InBounds(t)&&colors->Scale("RawQ").InBounds(q));

    // Only ghosted hits can be suppressed. Option to "ghost" dimmed 
    // hits to suppress them.
    if ( (drawopt->fRawDrawingOpt & evd::RawDrawingOptions::kSUPPRESS_DIMMED_HITS) && 
         DimFiltered(hit, unfiltered) ) *ghosted =  1;

    if (*ghosted | DimFiltered(hit, unfiltered ) ) {
      if(colors->Background(0) == kWhite) return 18; else return 12;
    }
    else {
      switch (drawopt->fColor) {
      case evd::RawDrawingOptions::kCOLOR_BY_CHARGE:
        co = q_color; 
        break;
      case evd::RawDrawingOptions::kCOLOR_BY_TIME:
        co = t_color; 
        break;
      default:
        co = 1;
      }
    }
    return co;
  }

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


  bool RawDataDrawer::DimFiltered(const rawdata::RawDigit* hit,
                                  const std::set<std::pair<uint32_t, int16_t>>& unfiltered)
  {
    static art::ServiceHandle<SliceNavigator> nav;
    if(!nav->DimHits() && !nav->DimOthers() ) return false;

    return unfiltered.count(std::make_pair(hit->Channel(), hit->TDC())) == 0;
  }

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

  ///
  /// Determine the scale factor to use when rendering hit cells
  ///
  double RawDataDrawer::ScaleFactor(double /*t*/, double q)
  {
    static art::ServiceHandle<evd::RawDrawingOptions> drawopt;
    
    double sf = drawopt->fScaleFactor;
    double q0 = 4096.0;
    if (drawopt->fRawDrawingOpt & evd::RawDrawingOptions::kSCALE_HITS_BY_CHARGE) {
      sf *= 3.0*sqrt(q/q0);
    }
    return sf;
  }
  
  ///
  /// Draw a stick figure to scale
  ///
  void RawDataDrawer::StickFigure(evdb::View3D* view,
                                  double x,
                                  double y,
                                  double z)
  {
    double h = 178.0; // Average height of a man
    double xlo= 9E9;
    double xhi=-9E9;
    double ylo= 9E9;
    double yhi=-9E9;
    double leg[8][2] = {
      {406,-618}, // 0
      {392,-576}, // 1
      {481,-509}, // 2
      {511,-403}, // 3
      {444,-385}, // 4
      {451,-495}, // 5
      {451,-618}, // 6
      {429,-617}  // 7
    };
    double arms[8][2] = {
      {537,-414}, // 0
      {535,-394}, // 1
      {552,-309}, // 2
      {550,-228}, // 3
      {462,-239}, // 4
      {401,-290}, // 5
      {453,-345}, // 6
      {447,-363}  // 7
    };
    double body[9][2] = {
      {479,-394}, // 0
      {504,-221}, // 1
      {493,-220}, // 2
      {478,-212}, // 3
      {473,-181}, // 4
      {491,-165}, // 5
      {521,-174}, // 6
      {522,-211}, // 7
      {504,-221}  // 8
    };
    int i;
    for (i=0; i<8; ++i) {
      if (leg[i][0]<xlo) xlo = leg[i][0];
      if (leg[i][0]>xhi) xhi = leg[i][0];
      if (leg[i][1]<ylo) ylo = leg[i][1];
      if (leg[i][1]>yhi) yhi = leg[i][1];
    }
    for (i=0; i<8; ++i) {
      if (arms[i][0]<xlo) xlo = arms[i][0];
      if (arms[i][0]>xhi) xhi = arms[i][0];
      if (arms[i][1]<ylo) ylo = arms[i][1];
      if (arms[i][1]>yhi) yhi = arms[i][1];
    }
    for (i=0; i<8; ++i) {
      if (body[i][0]<xlo) xlo = body[i][0];
      if (body[i][0]>xhi) xhi = body[i][0];
      if (body[i][1]<ylo) ylo = body[i][1];
      if (body[i][1]>yhi) yhi = body[i][1];
    }
    
    double sf = h/fabs(yhi-ylo);
    
    double z0 = z-sf*fabs(xhi-xlo)-2.0;
    double y0 = y+2.0;
    for (i=0; i<8; ++i) {
      leg[i][0] = z0+sf*(leg[i][0]-xlo);
      leg[i][1] = y0+sf*(leg[i][1]-ylo);
    }
    for (i=0; i<8; ++i) {
      arms[i][0] = z0+sf*(arms[i][0]-xlo);
      arms[i][1] = y0+sf*(arms[i][1]-ylo);
    }
    for (i=0; i<9; ++i) {
      body[i][0] = z0+sf*(body[i][0]-xlo);
      body[i][1] = y0+sf*(body[i][1]-ylo);
    }
    
    TPolyLine3D& l = view->AddPolyLine3D(8,kOrange-9,1,2);
    TPolyLine3D& a = view->AddPolyLine3D(8,kOrange-9,1,2);
    TPolyLine3D& b = view->AddPolyLine3D(9,kOrange-9,1,2);
    l.SetLineStyle(1);
    a.SetLineStyle(1);
    b.SetLineStyle(1);
    l.SetLineWidth(1);
    a.SetLineWidth(1);
    b.SetLineWidth(1);
    for (i=0; i<8; ++i) {
      l.SetPoint(i,x,leg[i][1],leg[i][0]);
    }
    for (i=0; i<8; ++i) {
      a.SetPoint(i,x,arms[i][1],arms[i][0]);
    }
    for (i=0; i<9; ++i) {
      b.SetPoint(i,x,body[i][1],body[i][0]);
    }
  }

} // end namespace evd
////////////////////////////////////////////////////////////////////////