packages/CaloReco/RawClusterBuilderTopo.h

0001 #ifndef CALORECO_RAWCLUSTERTOPO_H
0002 #define CALORECO_RAWCLUSTERTOPO_H
0003
0004 //===========================================================
0005 /// \file RawClusterBuilderTopo.h
0006 /// \brief 3-D topoClustering across calorimeter layers
0007 /// \author Dennis V. Perepelitsa
0008 //===========================================================
0009
0010 #include <fun4all/SubsysReco.h>
0011
0012 #include <map>
0013 #include <string>
0014 #include <utility>  // for pair
0015 #include <vector>
0016
0017 class PHCompositeNode;
0018 class RawClusterContainer;
0019 class RawTowerGeomContainer;
0020
0021 class RawClusterBuilderTopo : public SubsysReco
0022 {
0023  public:
0024   explicit RawClusterBuilderTopo(const std::string &name = "RawClusterBuilderTopo");
0025   ~RawClusterBuilderTopo() override {}
0026
0027   int InitRun(PHCompositeNode *topNode) override;
0028   int process_event(PHCompositeNode *topNode) override;
0029   int End(PHCompositeNode *topNode) override;
0030
0031   void set_nodename(const std::string &nodename)
0032   {
0033     ClusterNodeName = nodename;
0034   }
0035
0036   void set_noise(float noise_0 = 0.0025, float noise_1 = 0.006, float noise_2 = 0.03)
0037   {
0038     _noise_LAYER[0] = noise_0;
0039     _noise_LAYER[1] = noise_1;
0040     _noise_LAYER[2] = noise_2;
0041   }
0042
0043   void set_significance(float seed, float grow, float peri)
0044   {
0045     _sigma_seed = seed;
0046     _sigma_grow = grow;
0047     _sigma_peri = peri;
0048   }
0049
0050   void set_absE(bool allow)
0051   {
0052     _use_absE = allow;
0053   }
0054
0055   void allow_corner_neighbor(bool allow)
0056   {
0057     _allow_corner_neighbor = allow;
0058   }
0059
0060   void set_enable_HCal(bool enable_HCal)
0061   {
0062     _enable_HCal = enable_HCal;
0063   }
0064
0065   void set_enable_EMCal(bool enable_EMCal)
0066   {
0067     _enable_EMCal = enable_EMCal;
0068   }
0069
0070   void set_do_split(bool do_split)
0071   {
0072     _do_split = do_split;
0073   }
0074
0075   void set_minE_local_max(float minE_0 = 1, float minE_1 = 1, float minE_2 = 1)
0076   {
0077     _local_max_minE_LAYER[0] = minE_0;
0078     _local_max_minE_LAYER[1] = minE_1;
0079     _local_max_minE_LAYER[2] = minE_2;
0080   }
0081
0082   void set_R_shower(float R_shower)
0083   {
0084     _R_shower = R_shower;
0085   }
0086
0087   void set_use_only_good_towers(bool b)
0088   {
0089     _only_good_towers = b;
0090   }
0091   bool get_use_only_good_towers()
0092   {
0093     return _only_good_towers;
0094   }
0095
0096   void set_min_cluster_E_saved(float min_cluster_E)
0097   {
0098     _min_cluster_E = min_cluster_E;
0099   }
0100
0101  private:
0102   void CreateNodes(PHCompositeNode *topNode);
0103
0104   // geometric constants to express IHCal<->EMCal overlap in eta
0105   static int RawClusterBuilderTopo_constants_EMCal_eta_start_given_IHCal[];
0106
0107   static int RawClusterBuilderTopo_constants_EMCal_eta_end_given_IHCal[];
0108
0109   static int RawClusterBuilderTopo_constants_IHCal_eta_given_EMCal[];
0110
0111   // utility functions to express IHCal<->EMCal overlap in phi
0112   int get_first_matching_EMCal_phi_from_IHCal(int index_hcal_phi)
0113   {
0114     return (4 * index_hcal_phi + 5) % _EMCAL_NPHI;
0115   }
0116
0117   int get_matching_HCal_phi_from_EMCal(int index_emcal_phi)
0118   {
0119     return ((index_emcal_phi + 251) / 4) % _HCAL_NPHI;
0120   }
0121
0122   std::vector<int> get_adjacent_towers_by_ID(int ID);
0123
0124   static float calculate_dR(float, float, float, float);
0125
0126   void export_single_cluster(const std::vector<int> &);
0127
0128   void export_clusters(const std::vector<int> &, std::map<int, std::pair<int, int> >, unsigned int, const std::vector<float> &, const std::vector<float> &, const std::vector<float> &);
0129
0130   int get_ID(int ilayer, int ieta, int iphi)
0131   {
0132     if (ilayer < 2)
0133     {
0134       return ilayer * _HCAL_NETA * _HCAL_NPHI + ieta * _HCAL_NPHI + iphi;
0135     }
0136     return _EMCAL_NPHI * _EMCAL_NETA + ieta * _EMCAL_NPHI + iphi;
0137   }
0138
0139   int get_ilayer_from_ID(int ID)
0140   {
0141     if (ID < _EMCAL_NPHI * _EMCAL_NETA)
0142     {
0143       return ((int) (ID / (_HCAL_NETA * _HCAL_NPHI)));
0144     }
0145     else
0146     {
0147       return 2;
0148     }
0149   }
0150
0151   int get_ieta_from_ID(int ID)
0152   {
0153     if (ID < _EMCAL_NPHI * _EMCAL_NETA)
0154     {
0155       return ((int) ((ID % (_HCAL_NETA * _HCAL_NPHI)) / (_HCAL_NPHI)));
0156     }
0157     else
0158     {
0159       return ((int) ((ID - _EMCAL_NPHI * _EMCAL_NETA) / _EMCAL_NPHI));
0160     }
0161   }
0162
0163   int get_iphi_from_ID(int ID)
0164   {
0165     if (ID < _EMCAL_NPHI * _EMCAL_NETA)
0166     {
0167       return ((int) (ID % _HCAL_NPHI));
0168     }
0169     else
0170     {
0171       return ((int) ((ID - _EMCAL_NPHI * _EMCAL_NETA) % _EMCAL_NPHI));
0172     }
0173   }
0174
0175   int get_status_from_ID(int ID)
0176   {
0177     if (ID < _EMCAL_NPHI * _EMCAL_NETA)
0178     {
0179       return _TOWERMAP_STATUS_LAYER_ETA_PHI[get_ilayer_from_ID(ID)][get_ieta_from_ID(ID)][get_iphi_from_ID(ID)];
0180     }
0181     return _EMTOWERMAP_STATUS_ETA_PHI[get_ieta_from_ID(ID)][get_iphi_from_ID(ID)];
0182   }
0183
0184   float get_E_from_ID(int ID)
0185   {
0186     if (ID < _EMCAL_NPHI * _EMCAL_NETA)
0187     {
0188       return _TOWERMAP_E_LAYER_ETA_PHI[get_ilayer_from_ID(ID)][get_ieta_from_ID(ID)][get_iphi_from_ID(ID)];
0189     }
0190     return _EMTOWERMAP_E_ETA_PHI[get_ieta_from_ID(ID)][get_iphi_from_ID(ID)];
0191   }
0192
0193   void set_status_by_ID(int ID, int status)
0194   {
0195     if (ID < _EMCAL_NPHI * _EMCAL_NETA)
0196     {
0197       _TOWERMAP_STATUS_LAYER_ETA_PHI[get_ilayer_from_ID(ID)][get_ieta_from_ID(ID)][get_iphi_from_ID(ID)] = status;
0198     }
0199     else
0200     {
0201       _EMTOWERMAP_STATUS_ETA_PHI[get_ieta_from_ID(ID)][get_iphi_from_ID(ID)] = status;
0202     }
0203   }
0204
0205   RawClusterContainer *_clusters {nullptr};
0206
0207   RawTowerGeomContainer *_geom_containers[3]{};
0208
0209   // geometric parameters defined at runtime
0210   int _EMCAL_NETA {-1};
0211   int _EMCAL_NPHI{-1};
0212
0213   int _HCAL_NETA{-1};
0214   int _HCAL_NPHI{-1};
0215
0216   float _noise_LAYER[3]{};
0217
0218   float _sigma_seed {4.0};
0219   float _sigma_grow {2.0};
0220   float _sigma_peri {0.0};
0221   float _local_max_minE_LAYER[3]{};
0222   float _R_shower {0.025};
0223   float _min_cluster_E{0.0};
0224
0225   bool _allow_corner_neighbor {true};
0226   bool _use_absE {true};
0227
0228   bool _enable_HCal {true};
0229   bool _enable_EMCal {true};
0230
0231   bool _do_split {true};
0232   bool _only_good_towers {true};
0233
0234   std::vector<std::vector<std::vector<float> > > _TOWERMAP_E_LAYER_ETA_PHI;
0235   std::vector<std::vector<std::vector<int> > > _TOWERMAP_KEY_LAYER_ETA_PHI;
0236   std::vector<std::vector<std::vector<int> > > _TOWERMAP_STATUS_LAYER_ETA_PHI;
0237
0238   std::vector<std::vector<float> > _EMTOWERMAP_E_ETA_PHI;
0239   std::vector<std::vector<int> > _EMTOWERMAP_KEY_ETA_PHI;
0240   std::vector<std::vector<int> > _EMTOWERMAP_STATUS_ETA_PHI;
0241
0242   std::string ClusterNodeName {"TOPOCLUSTER_HCAL"};
0243 };
0244
0245 #endif