From ba785f450c91484cf0e3cbfcfb1aab356aa01614 Mon Sep 17 00:00:00 2001 From: Marco Mirazita Date: Fri, 26 Jun 2026 03:08:57 -0400 Subject: [PATCH] Rich clustering update --- .../detectors/RICHBankBuilder.java | 8 +- .../java/org/jlab/rec/rich/RICHCluster.java | 55 ++++++-- .../main/java/org/jlab/rec/rich/RICHHit.java | 2 + .../jlab/rec/rich/RICHPMTReconstruction.java | 128 +++++++++++++----- 4 files changed, 143 insertions(+), 50 deletions(-) diff --git a/reconstruction/calib/src/main/java/org/jlab/calibration/detectors/RICHBankBuilder.java b/reconstruction/calib/src/main/java/org/jlab/calibration/detectors/RICHBankBuilder.java index 4309ff89ef..3a4228f389 100644 --- a/reconstruction/calib/src/main/java/org/jlab/calibration/detectors/RICHBankBuilder.java +++ b/reconstruction/calib/src/main/java/org/jlab/calibration/detectors/RICHBankBuilder.java @@ -20,9 +20,10 @@ public RICHBankBuilder() { } @Override public boolean isGoodEvent(DataEvent event) { + + DataBank part = event.getBank("REC::Particle"); - DataBank part = event.getBank("REC::Particle"); - if(part.rows()<1 || + if(part.rows()<1 || part.getInt("pid", 0)!=11 || ((int) (Math.abs(part.getShort("status", 0))/1000))!=2) return false; @@ -37,7 +38,8 @@ public DataBank buildCalibBank(DataEvent event) { DataBank hits = event.getBank("RICH::Hit"); DataBank clus = event.getBank("RICH::Cluster"); DataBank phos = event.getBank("RICH::Photon"); - + + List goodPhotons = new ArrayList<>(); for(int i=0; i { +public class RICHCluster extends ArrayList implements Comparable { /** * A cluster in the RICH consists of an array of anodes in one PMT @@ -339,14 +339,31 @@ public void merge(RICHCluster clu) { } + // ---------------- + public boolean overlaps(RICHCluster clu) { + // ---------------- + + if (this.get(0).get_sector() != clu.get(0).get_sector() ) return false; + if (this.get(0).get_pmt() != clu.get(0).get_pmt() ) return false; + + for(int i=0; i comp "+this.get_channel()+" "+this.get_charge()+" "+ocluster.get_channel()+" "+ocluster.get_charge()); - if(this.get_charge() == ocluster.get_charge())return 0; - if(this.get_charge() > ocluster.get_charge()){ - return 1; - }else{ + //if(this.get_charge() == ocluster.get_charge())return 0; + //if(this.get_charge() > ocluster.get_charge()){ + // return 1; + //}else{ + // return -1; + //} + + if(this.get(0).get_sector()!=ocluster.get(0).get_sector()) + return this.get(0).get_sector()>ocluster.get(0).get_sector() ? 1 : -1; + + if(this.get(0).get_pmt()!=ocluster.get(0).get_pmt()) + return this.get(0).get_pmt()>ocluster.get(0).get_pmt() ? 1 : -1; + + if(this.get_size() == ocluster.get_size()) + return 0; + else if(this.get_size() > ocluster.get_size()) return -1; - } + else + return 1; + + } // ---------------- public void showCluster() { // ---------------- - System.out.format("Cluster ID %3d PMT %4d Siz %4d Ch %7.1f T %7.1f raw %7.1f wT %7.1f glxy %4d %4d XYZ %7.2f %7.2f %7.2f wXYZ %7.2f %7.2f %7.2f \n", + System.out.format("Cluster ID %3d sector %1d PMT %4d Siz %4d Ch %7.1f T %7.1f raw %7.1f wT %7.1f glxy %4d %4d XYZ %7.2f %7.2f %7.2f wXYZ %7.2f %7.2f %7.2f \n", this.clusid, + this.get(0).get_sector(), this.get(0).get_pmt(), this.get_size(), this.get_charge(), @@ -383,9 +418,9 @@ public void showCluster() { this.get_x(), this.get_y(), this.get_z(), this.get_wx(), this.get_wy(), this.get_wz()); for(int j = 0; j< this.size(); j++) { - System.out.format(" --> hit # %3d ID %3d idxy %3d %3d dur %4d \n", - j, this.get(j).get_id(), - this.get(j).get_idx(), this.get(j).get_idy(), this.get(j).get_duration()); + System.out.format(" --> hit # %3d ID %3d anode %2d idxy %3d %3d dur %4d clu %2d\n", + j, this.get(j).get_id(), this.get(j).get_anode(), + this.get(j).get_idx(), this.get(j).get_idy(), this.get(j).get_duration(), this.get(j).get_cluster()); } } diff --git a/reconstruction/rich/src/main/java/org/jlab/rec/rich/RICHHit.java b/reconstruction/rich/src/main/java/org/jlab/rec/rich/RICHHit.java index 67ac427280..5de2f9ccd2 100644 --- a/reconstruction/rich/src/main/java/org/jlab/rec/rich/RICHHit.java +++ b/reconstruction/rich/src/main/java/org/jlab/rec/rich/RICHHit.java @@ -366,6 +366,8 @@ public boolean passHitSelection(RICHHit hit) { } + // ---------------- + @Override // ---------------- public int compareTo(RICHHit ohit) { // ---------------- diff --git a/reconstruction/rich/src/main/java/org/jlab/rec/rich/RICHPMTReconstruction.java b/reconstruction/rich/src/main/java/org/jlab/rec/rich/RICHPMTReconstruction.java index 31ad6280ae..ecbd39e728 100644 --- a/reconstruction/rich/src/main/java/org/jlab/rec/rich/RICHPMTReconstruction.java +++ b/reconstruction/rich/src/main/java/org/jlab/rec/rich/RICHPMTReconstruction.java @@ -270,10 +270,12 @@ public ArrayList findClusters(ArrayList hits) { // ---------------- int debugMode = 0; + + ArrayList allclusters = new ArrayList(); if(debugMode>=2) { System.out.println("--------------------\n"); - System.out.println("Building allclusters\n"); + System.out.println("Building allclusters for event number " + richevent.get_EventID() + "\n"); System.out.println("--------------------\n"); } @@ -301,12 +303,16 @@ public ArrayList findClusters(ArrayList hits) { } } - if(debugMode>=2){ + // Ordering the cluster + Collections.sort(allclusters); + + + if(debugMode>=2){ System.out.println("List of all Clusters"); for(int i=0; i selectGoodClusters(ArrayList allclust // ---------------- int debugMode = 0; + + if(debugMode>=1) { + System.out.println("--------------------\n"); + System.out.println("Selecting good clusters"); + System.out.println("--------------------\n"); + } + ArrayList clusters = new ArrayList(); - if(debugMode>=2) System.out.println("\nSelecting good clusters"); + ArrayList goodclusters = new ArrayList(); + int nclu = 0; for(int i=0; i=1)System.out.format(" merge clu %d %d xyz %7.2f %7.2f %7.2f to %d %d \n",i,goodclu.get_id(),goodclu.get_x(), goodclu.get_y(), goodclu.get_z(), j,clusters.get(j).get_id()); - } - } - if(merge==0){ - nclu++; - goodclu.set_id(nclu); - clusters.add(goodclu); - } - }else{ - // cancel hit to cluster link - RICHCluster badclu = allclusters.get(i); - for(int j = 0; j< badclu.size(); j++) { - badclu.get(j).set_cluster(0); - } - } - } + RICHCluster clu = allclusters.get(i); - if(debugMode>=1){ + if (debugMode>=1) { + System.out.println(" ---> Looking at cluster " + i + " id=" + clu.get_id() + " Sector=" + clu.get(0).get_sector() + " PMT=" + clu.get(0).get_pmt() + " size=" + clu.size()); + for (int k=0; k=1) { + System.out.println(" merging cluster " + clu.get_id() + " to cluster " + clusters.get(j).get_id() + " newsize=" + clusters.get(j).size() ); + } + } + } + + + if(merge==0){ + nclu++; + clu.set_id(nclu); + clusters.add(clu); + + if(debugMode>=1) { + System.out.println(" Creating new cluster " + clusters.get(nclu-1).get_id() + " size=" + clusters.get(nclu-1).size() ); + } + + } + + } + + /*if(debugMode>=1){ System.out.format("-------------------------\n"); - System.out.format("List of selected Clusters %4d \n",clusters.size()); + System.out.format("List of merged Clusters %4d \n",clusters.size()); System.out.format("-------------------------\n"); for(int i=0; i=1) { + System.out.println(" ---> Selecting good cluster " + i + " id=" + clu.get_id() + " Sector=" + clu.get(0).get_sector() + " PMT=" + clu.get(0).get_pmt() + " size=" + clu.size()); + } - /*Collections.sort(clusters); + ngoodclu++; + clu.set_id(ngoodclu); + for (int j=0; j=1){ - System.out.println("List of sorted Clusters"); - for(int i=0; i hits, ArrayList allcluste if(hiti.get_cluster()!=0) continue; // this hit is not yet associated with a cluster if(debugMode==1)System.out.println("Hit pair "+ih+" "+hiti.get_id()+" "+hiti.get_pmt()+" "+hiti.get_channel()+" "+hiti.get_duration()+" "+hiti.get_cluster()+" | " +jh+" "+hitj.get_id()+" "+hitj.get_pmt()+" "+hitj.get_channel()+" "+hitj.get_duration()+" "+hitj.get_cluster()); - if(hiti.get_pmt()==hitj.get_pmt() && hitj.get_duration()*100 < hiti.get_duration()*richpar.GOODHIT_FRAC){ + if(hiti.get_sector()==hitj.get_sector() && hiti.get_pmt()==hitj.get_pmt() && hitj.get_duration()*100 < hiti.get_duration()*richpar.GOODHIT_FRAC){ for(int k=-1; k<=1; k+=2 ) { if(hiti.get_channel() == (k+hitj.get_channel())) {hitj.set_xtalk(1000+hiti.get_id()+1); if(debugMode==1)System.out.println(" E Xtalk "+hitj.get_xtalk());} }