From 631bc0e88b518b37cb5e15bc730258cdee838d55 Mon Sep 17 00:00:00 2001 From: raffaelladevita Date: Fri, 21 Nov 2025 19:35:04 -0500 Subject: [PATCH 01/15] switch phi binning to R1 phi and improved RG-L vertex fit --- .../java/org/clas/dc/alignment/Histo.java | 51 +++++++++++-------- 1 file changed, 30 insertions(+), 21 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 7516f33..32aadf2 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -353,6 +353,12 @@ public double phiSector() { public double thetaDC() { return this.traj.theta(); } + + public double phiDC() { + Vector3D dir = new Vector3D(this.traj); + dir.rotateZ(-Math.PI/3*(sector-1)); + return dir.phi(); + } } private int getElectronIndex(Event event) { @@ -473,7 +479,7 @@ private void processEvent(Event event, Event shifted) { electron.vector().rotateZ(Math.toRadians(-60*(electron.sector()-1))); double theta = Math.toDegrees(electron.thetaDC()); - double phi = Math.toDegrees(electron.phi()); + double phi = Math.toDegrees(electron.phiDC()); double vz = electron.vz(); int sector = electron.sector(); @@ -1638,19 +1644,22 @@ public static void fitRGLVertex(H1F histo) { int nbin = histo.getData().length; double dx = histo.getDataX(1)-histo.getDataX(0); //find downstream window - int ibin0 = Histo.getMaximumBinBetween(histo, histo.getDataX(0), (Constants.TARGETPOS-Constants.TARGETLENGTH/2)); + int ibin1 = Histo.getMaximumBinBetween(histo, histo.getDataX(0), Constants.TARGETCENTER); + int ibin2 = Histo.getMaximumBinBetween(histo, Constants.TARGETCENTER, histo.getDataX(nbin-1)); //check if the found maximum is the first or second peak, ibin is tentative upstream window - int ibin1 = Math.max(0, ibin0 - (int)(Constants.TARGETLENGTH/dx)); - int ibin2 = Math.min(nbin-1, ibin0 + (int)(Constants.TARGETLENGTH/dx)); - if(histo.getBinContent(ibin1) Date: Mon, 1 Dec 2025 19:06:20 -0500 Subject: [PATCH 02/15] bug fix in region summary graphs --- .../src/main/java/org/clas/dc/alignment/Alignment.java | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java index 16c3b78..6250ae1 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java @@ -865,7 +865,7 @@ private DataGroup getRegionHistograms(String parameter, Table alignment, int ico residuals.addDataSet(hi_res, ir*Constants.NSECTOR + is); for(int it=1; it Date: Fri, 5 Dec 2025 15:09:47 -0500 Subject: [PATCH 03/15] lowering logging level to INFO, correct usage printout for vertex fits, saved vertex fit parameters for recent datasets, rotated wire vs. time plots --- .../java/org/clas/dc/alignment/Alignment.java | 23 +++++++++++++------ .../java/org/clas/dc/alignment/Constants.java | 5 +++- .../java/org/clas/dc/alignment/Histo.java | 10 ++++---- 3 files changed, 26 insertions(+), 12 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java index 6250ae1..f387e8e 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java @@ -201,6 +201,12 @@ else if(Math.abs(vertex)==6) Constants.initTargetPars(Constants.RGCSUMMER2022); else if(Math.abs(vertex)==7) Constants.initTargetPars(Constants.RGDFALL2023); + else if(Math.abs(vertex)==8) + Constants.initTargetPars(Constants.RGASPRING2018); + else if(Math.abs(vertex)==9) + Constants.initTargetPars(Constants.RGKSPRING2024); + else if(Math.abs(vertex)==10) + Constants.initTargetPars(Constants.RGESPRING2024); } if(vertex<0) vertex=0; } @@ -1169,7 +1175,7 @@ public void saveHistos(String fileName) { public static void main(String[] args){ Alignment align = new Alignment(); - align.initLogger(Level.CONFIG); + align.initLogger(Level.INFO); String[] inputs = align.getInputs(); @@ -1189,9 +1195,10 @@ public static void main(String[] args){ parser.getOptionParser("-process").addOption("-time" , "0", "make time residual histograms (1=true, 0=false)"); parser.getOptionParser("-process").addOption("-residuals", "2", "fit residuals with double gaussian (2), single gaussian (1), or use mean (0)"); parser.getOptionParser("-process").addOption("-vertfit" , "5", "fit vertex plots with:\n" + + "\t\t- RG-L layout (11), \n" + "\t\t- RG-E layout (10), dual target, \n" + - "\t\t- RG-A Spring18 layout (9)\n" + - "\t\t- RG-K layout (8), new cryotarget, \n" + + "\t\t- RG-K layout (9), new cryotarget, \n" + + "\t\t- RG-A Spring18 layout (8)\n" + "\t\t- RG-D layout (7), new cryotarget, \n" + "\t\t- RG-C layout (6),\n" + "\t\t- 4 gaussians (5),\n" + @@ -1228,9 +1235,10 @@ public static void main(String[] args){ parser.getOptionParser("-analyze").addOption("-shifts" , "0", "use event-by-event subtraction for unit shifts (1=on, 0=off)"); parser.getOptionParser("-analyze").addOption("-residuals", "2", "fit residuals (2), use mean (1), or use existing fit available (0)"); parser.getOptionParser("-analyze").addOption("-vertfit" , "5", "fit vertex plots with:\n" + + "\t\t- RG-L layout (11), \n" + "\t\t- RG-E layout (10), dual target, \n" + - "\t\t- RG-A Spring18 layout (9)\n" + - "\t\t- RG-K layout (8), new cryotarget, \n" + + "\t\t- RG-K layout (9), new cryotarget, \n" + + "\t\t- RG-A Spring18 layout (8)\n" + "\t\t- RG-D layout (7), new cryotarget,\n" + "\t\t- RG-C layout (6),\n" + "\t\t- 4 gaussians (5),\n" + @@ -1262,9 +1270,10 @@ public static void main(String[] args){ parser.getOptionParser("-fit").addOption("-shifts" , "0", "use event-by-event subtraction for unit shifts (1=on, 0=off)"); parser.getOptionParser("-fit").addOption("-sector" , "1", "sector-dependent derivatives (1) or average (0)"); parser.getOptionParser("-fit").addOption("-vertfit" , "5", "fit vertex plots with:\n" + + "\t\t- RG-L layout (11), \n" + "\t\t- RG-E layout (10), dual target, \n" + - "\t\t- RG-A Spring18 layout (9)\n" + - "\t\t- RG-K layout (8), new cryotarget, \n" + + "\t\t- RG-K layout (9), new cryotarget, \n" + + "\t\t- RG-A Spring18 layout (8)\n" + "\t\t- RG-D layout (7), new cryotarget, \n" + "\t\t- RG-C layout (6),\n" + "\t\t- 4 gaussians (5),\n" + diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java index e3759fd..52b9668 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java @@ -89,12 +89,15 @@ public class Constants { // - distance between cell exit window and indulation foil, // - distance beetween the scattering chamber exit window and the target center public static final double[] DEFAULT = {-0.5, 5.0, 6.8, 27.3}; + public static final double[] RGASPRING2018 = { 0.0, 5.0, 2.8, 28.6}; public static final double[] RGAFALL2018 = {-0.5, 5.0, 2.8, 28.4}; public static final double[] RGBSPRING2019 = {-0.5, 5.0, 6.8, 28.4}; public static final double[] RGFSUMMER2020 = {-32, 5.0, 6.8, 27.3}; public static final double[] RGMFALL2021 = {-0.5, 5.0, 6.8, 27.3}; public static final double[] RGCSUMMER2022 = {-1.4,5.25, 8.3, 14.3}; - public static final double[] RGDFALL2023 = {-2.5, 5.0, 3.0, 27.1}; + public static final double[] RGDFALL2023 = {-3.0, 5.0, 3.0, 27.1}; + public static final double[] RGKSPRING2024 = {-0.5, 5.0, 0.0, 26.5}; // warm + public static final double[] RGESPRING2024 = {-5.5, 2.2, 2.0, 25.9}; // target parameters used for vertex fit initialization public static double TARGETPOS = DEFAULT[0]; diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 32aadf2..e8f9824 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -194,9 +194,9 @@ private void createHistos(String optStats) { this.wires[is] = new DataGroup(nSector, nSLayer); for(int il=0; il Date: Mon, 26 Jan 2026 16:32:49 -0500 Subject: [PATCH 04/15] fixed vertex range for non-nominal histos --- dc/java/dc-alignment/pom.xml | 2 +- .../java/org/clas/dc/alignment/Histo.java | 2 + dc/utilities/elastic.groovy | 843 ++++++++++++++++++ 3 files changed, 846 insertions(+), 1 deletion(-) create mode 100644 dc/utilities/elastic.groovy diff --git a/dc/java/dc-alignment/pom.xml b/dc/java/dc-alignment/pom.xml index 8c9ebaa..3807c0a 100644 --- a/dc/java/dc-alignment/pom.xml +++ b/dc/java/dc-alignment/pom.xml @@ -57,7 +57,7 @@ jhep-maven - https://clasweb.jlab.org/jhep/maven + https://clasweb.jlab.org/.jhep/maven diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index e8f9824..536b72c 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -77,6 +77,8 @@ public Histo(String name, List files, Bin[] thetabins, Bin[] phibins, do this.name = name; this.thetaBins = thetabins; this.phiBins = phibins; + this.minVtx = vertexrange[0]; + this.maxVtx = vertexrange[1]; this.nominalFiles = files; this.createHistos(optstats); } diff --git a/dc/utilities/elastic.groovy b/dc/utilities/elastic.groovy new file mode 100644 index 0000000..3595fb7 --- /dev/null +++ b/dc/utilities/elastic.groovy @@ -0,0 +1,843 @@ +package org.clas.dc.alignment; + +import java.io.File; +import java.util.LinkedHashMap; +import java.util.List; +import java.util.Map; +import java.util.logging.Level; +import java.util.logging.Logger; +import javax.swing.JFrame; +import org.jlab.clas.pdg.PhysicsConstants; +import org.jlab.clas.physics.LorentzVector; +import org.jlab.clas.physics.Particle; +import org.jlab.clas.physics.PhysicsEvent; +import org.jlab.groot.base.GStyle; +import org.jlab.groot.data.DataLine; +import org.jlab.groot.data.H1F; +import org.jlab.groot.data.H2F; +import org.jlab.groot.data.IDataSet; +import org.jlab.groot.data.TDirectory; +import org.jlab.groot.fitter.DataFitter; +import org.jlab.groot.graphics.EmbeddedCanvasTabbed; +import org.jlab.groot.graphics.EmbeddedPad; +import org.jlab.groot.graphics.IDataSetPlotter; +import org.jlab.groot.group.DataGroup; +import org.jlab.groot.math.F1D; +import org.jlab.io.base.DataBank; +import org.jlab.io.base.DataEvent; +import org.jlab.io.hipo.HipoDataSource; +import org.jlab.jnp.utils.benchmark.ProgressPrintout; +import org.jlab.jnp.utils.options.OptionStore; + +/** + * + * @author devita + */ +public class Kinematics { + + private double ebeam = 10.6; + private double targetPos = -3.5; + private double targetLength = 5; +// private double scWindow = 17.5+targetPos; + private double scWindow = 28.4+targetPos; + private String fontName = "Arial"; + private File lundfile = null; + + private Map histos = new LinkedHashMap<>(); + + + public Kinematics(double energy, String vertex) { + this.initGraphics(); + this.initVertexPar(vertex); + this.setEbeam(energy); + } + + private void initVertexPar(String pars) { + if(!pars.isEmpty()) { + double[] parValues = new double[pars.split(":").length]; + for(int i=0; i0) targetPos = parValues[0]; + if(parValues.length>1) targetLength = parValues[1]; + if(parValues.length>2) scWindow = parValues[2]+targetPos; + } + } + + public double getEbeam() { + return ebeam; + } + + public void setEbeam(double ebeam) { + this.ebeam = ebeam; + } + + public double getMaxQ2() { + double maxQ2 = 3; + if(ebeam>6) maxQ2 = 5; + return maxQ2; + } + + public double getMaxW() { + double maxW = ebeam*0.4+0.6; + if(ebeam>6.5) maxW = 4.5; + return maxW; + } + + private void initGraphics() { + GStyle.getAxisAttributesX().setTitleFontSize(18); + GStyle.getAxisAttributesX().setLabelFontSize(14); + GStyle.getAxisAttributesY().setTitleFontSize(18); + GStyle.getAxisAttributesY().setLabelFontSize(14); + GStyle.getAxisAttributesZ().setLabelFontSize(14); + GStyle.getAxisAttributesX().setLabelFontName(this.fontName); + GStyle.getAxisAttributesY().setLabelFontName(this.fontName); + GStyle.getAxisAttributesZ().setLabelFontName(this.fontName); + GStyle.getAxisAttributesX().setTitleFontName(this.fontName); + GStyle.getAxisAttributesY().setTitleFontName(this.fontName); + GStyle.getAxisAttributesZ().setTitleFontName(this.fontName); + GStyle.setGraphicsFrameLineWidth(1); + GStyle.getH1FAttributes().setLineWidth(2); + } + + public void createHistos() { + // General + H2F hi_q2w = new H2F("hi_q2w","hi_q2w",100, 0.6, this.getMaxW(), 100, 0.0, this.getMaxQ2()); + hi_q2w.setTitleX("W (GeV)"); + hi_q2w.setTitleY("Q2 (GeV2)"); + H1F hi_w = new H1F("hi_w","hi_w",500, 0.6, this.getMaxW()); + hi_w.setTitleX("W (GeV)"); + hi_w.setTitleY("Counts"); + H2F hi_w_phi = new H2F("hi_w_phi","hi_w_phi",100, -180.0, 180.0, 250, 0.6, this.getMaxW()); + hi_w_phi.setTitleX("#phi (deg)"); + hi_w_phi.setTitleY("W (GeV)"); + H2F hi_z_phi = new H2F("hi_z_phi","hi_z_phi",100, -180.0, 180.0, 100, 20, 30); + hi_z_phi.setTitleX("#phi (deg)"); + hi_z_phi.setTitleY("z (cm)"); + H2F hi_w_theta = new H2F("hi_w_theta","hi_w_theta",100, 5.0, 15.0, 250, 0.6, this.getMaxW()); + hi_w_theta.setTitleX("#theta (deg)"); + hi_w_theta.setTitleY("W (GeV)"); + H2F hi_el = new H2F("hi_el","hi_el",100, 0.5, ebeam+0.5, 100, 0.0, 35.0); + hi_el.setTitleX("p (GeV)"); + hi_el.setTitleY("#theta (deg)"); + hi_el.setTitle("Electron"); + F1D f1_el = new F1D("f1_el", "2*(180/3.14)*atan(sqrt(0.93832*([e0]-x)/2/[e0]/x))", ebeam*0.75, ebeam*0.99); + f1_el.setParameter(0, ebeam); + DataGroup dg_general = new DataGroup(3,2); + dg_general.addDataSet(hi_q2w, 0); + dg_general.addDataSet(hi_w_theta, 1); + dg_general.addDataSet(hi_el, 2); + dg_general.addDataSet(f1_el, 2); + dg_general.addDataSet(hi_w, 3); + dg_general.addDataSet(hi_w_phi, 4); + dg_general.addDataSet(hi_z_phi, 5); + histos.put("General", dg_general); + // Elastic + H1F hi_ela = new H1F("hi_w","hi_w",500, 0.6, 1.3); + hi_ela.setTitleX("W (GeV)"); + hi_ela.setTitleY("Counts"); + H2F hi_ela_phi = new H2F("hi_w_phi","hi_w_phi",100, -180.0, 180.0, 250, 0.6, 1.3); + hi_ela_phi.setTitleX("#phi (deg)"); + hi_ela_phi.setTitleY("W (GeV)"); + H2F hi_ela_theta = new H2F("hi_w_theta","hi_w_theta",100, 5.0, 15.0, 250, 0.6, 1.3); + hi_ela_theta.setTitleX("#theta (deg)"); + hi_ela_theta.setTitleY("W (GeV)"); + H2F hi_ela_p = new H2F("hi_w_p","hi_w_p",100, ebeam*0.8, ebeam, 100, 0.6, 1.3); + hi_ela_p.setTitleX("p (GeV)"); + hi_ela_p.setTitleY("W (GeV)"); + hi_ela_p.setTitle("Electron"); + DataGroup dg_elastic = new DataGroup(2,2); + dg_elastic.addDataSet(hi_ela, 0); + dg_elastic.addDataSet(hi_ela_p, 1); + dg_elastic.addDataSet(hi_ela_theta, 2); + dg_elastic.addDataSet(hi_ela_phi, 3); + histos.put("Elastic", dg_elastic); + // vertex + String[] names= ["pi+", "pi-", "el"]; + DataGroup dg_vertex = new DataGroup(3, 4); + for(int i=0; i1) rmin = -0.5; + H1F hi_mmass = new H1F("mxt" + i, "", 200, rmin, rmax); + hi_mmass.setTitleX(a2pimx[i-1]); + hi_mmass.setTitleY("Counts"); + // hi_mmass.setLineColor(col); + H1F hi_imass = new H1F("mit" + i, "", 200, 0.0, 3.0); + hi_imass.setTitleX("M(#pi#pi) (GeV)"); + hi_imass.setTitleY("Counts"); + // hi_imass.setLineColor(col); + dg_2pi.addDataSet(hi_mmass, i-1); + dg_2pi.addDataSet(hi_imass, i+2); + } + histos.put("2pi", dg_2pi); + // 1 pi + DataGroup dg_1pi = new DataGroup(2,2); + String[] atitle = ["M_X(ep#rarrow e'^#pi+X) (GeV)","M_^X2(ep#rarrow e'pX) (Ge^V2)"]; + for(int i=1; i<=2; i++) { + double rmin = 0.5; + double rmax = 2; + if(i>1) rmin = -0.5; + H2F hi_mw = new H2F("W" + i, "", 100, -180, 180, 100, rmin, rmax); + hi_mw.setTitleX("W (GeV)"); + hi_mw.setTitleY(atitle[i-1]); + // hi_w.setLineColor(col); + H1F hi_mmass = new H1F("mxt" + i, "", 400, rmin, rmax); + hi_mmass.setTitleX(atitle[i-1]); + hi_mmass.setTitleY("Counts"); + // hi_mmass.setLineColor(col); + dg_1pi.addDataSet(hi_mw, i-1); + dg_1pi.addDataSet(hi_mmass, i+1); + } + histos.put("1pi", dg_1pi); + // Proton + DataGroup dg_proton = new DataGroup(4,2); + H2F hi_pr_p_dp = new H2F("hi_p_dp", "", 100, 1, 4, 100, ebeam*0.75, ebeam*1.2); + hi_pr_p_dp.setTitleX("p (GeV)"); hi_pr_p_dp.setTitleY("Ebeam (GeV)"); + H2F hi_pr = new H2F("hi_pr","hi_pr",100, 0, ebeam/2-0.1, 100, 20.0, 80.0); + hi_pr.setTitleX("p (GeV)"); + hi_pr.setTitleY("#theta (deg)"); + hi_pr.setTitle("Proton"); + F1D f1_pr = new F1D("f1_pr", "(180/3.14)*acos(([e0]*[e0]+x*x-pow(([e0]+0.93832-sqrt(x*x+0.9382*0.9382)),2))/2/[e0]/x)", ebeam*0.08, ebeam*0.4); + f1_pr.setParameter(0, ebeam); + H2F hi_phi = new H2F("hi_phi", "hi_phi", 200, -180, 180, 200, -180, 180); + hi_phi.setTitleX("El #phi (deg)"); + hi_phi.setTitleY("Pr #phi (deg)"); + // hi_dphi.setOptStat("1110"); + H2F hi_dphi = new H2F("hi_dphi", "hi_dphi", 200, -180, 180, 200, -10, 10); + hi_dphi.setTitleX("Pr #phi (deg)"); + hi_dphi.setTitleY("#Delta#phi (deg)"); + H2F hi_dpr = new H2F("hi_dpr", "hi_dpr", 200, -180, 180, 200, -10, 10); + hi_dpr.setTitleX("Pr #phi (deg)"); + hi_dpr.setTitleY("#Delta p (GeV)"); + H2F hi_pr_phi_dphi = new H2F("hi_phi_dphi", "", 100, -180.0, 180.0, 100, 170.0, 190.0); + hi_pr_phi_dphi.setTitleX("#phi (deg)"); hi_pr_phi_dphi.setTitleY("#Delta#phi (deg)"); + H2F hi_pr_phi_dz = new H2F("hi_phi_dz", "", 100, -180.0, 180.0, 100, -5.0, 5.0); + hi_pr_phi_dz.setTitleX("#phi (deg)"); hi_pr_phi_dz.setTitleY("#Deltaz (cm)"); + H1F hi_pr_dp = new H1F("hi_dp", "Ebeam (GeV)", "Counts", 100, ebeam*0.75, ebeam*1.2); + H1F hi_pr_dtheta = new H1F("hi_dtheta", "#Delta#theta (deg)", "Counts", 100, -15.0, 15.0); + H1F hi_pr_dphi = new H1F("hi_dphi", "#Delta#phi (deg)", "Counts", 100, 170.0, 190.0); + H1F hi_pr_dz = new H1F("hi_dz", "#Deltaz (cm)", "Counts", 100, -5.0, 5.0); + dg_proton.addDataSet(hi_pr_p_dp, 0); + dg_proton.addDataSet(hi_pr, 1); + dg_proton.addDataSet(f1_pr, 1); + dg_proton.addDataSet(hi_pr_phi_dphi, 2); + dg_proton.addDataSet(hi_pr_phi_dz, 3); + dg_proton.addDataSet(hi_pr_dp, 4); + dg_proton.addDataSet(hi_pr_dtheta, 5); + dg_proton.addDataSet(hi_pr_dphi, 6); + dg_proton.addDataSet(hi_pr_dz, 7); + histos.put("Proton", dg_proton); + // Phi + DataGroup dg_phi = new DataGroup(2,3); + for(int sector=1; sector <= 6; sector++) { + H1F hi_dphi_sec = new H1F("hi_dphi_" + sector, "hi_dphi_" + sector, 100, 160.0, 200.0); + hi_dphi_sec.setTitleX("#Delta#phi (deg)"); + hi_dphi_sec.setTitleY("Counts"); + hi_dphi_sec.setTitle("Sector " + sector); + dg_phi.addDataSet(hi_dphi_sec, sector-1); + } + histos.put("Phi", dg_phi); + // Beam + DataGroup dg_beam = new DataGroup(2,3); + for(int sector=1; sector <= 6; sector++) { + H1F hi_beam_sec = new H1F("hi_beam_" + sector, "hi_beam_" + sector, 100, ebeam*0.75, ebeam*1.2); + hi_beam_sec.setTitleX("Beam Energy (GeV)"); + hi_beam_sec.setTitleY("Counts"); + hi_beam_sec.setTitle("Sector " + sector); + dg_beam.addDataSet(hi_beam_sec, sector-1); + } + histos.put("Beam", dg_beam); + } + + public void processEvent(DataEvent event) { + + DataBank recEvent = null; + DataBank recPart = null; + DataBank recScint = null; + DataBank recCal = null; + DataBank recTracks = null; + DataBank recTraj = null; + + if(event.hasBank("REC::Event")) recEvent = event.getBank("REC::Event"); + if(event.hasBank("REC::Particle")) recPart = event.getBank("REC::Particle"); + if(event.hasBank("REC::Scintillator")) recScint = event.getBank("REC::Scintillator"); + if(event.hasBank("REC::Calorimeter")) recCal = event.getBank("REC::Calorimeter"); + if(event.hasBank("REC::Track")) recTracks = event.getBank("REC::Track"); + if(event.hasBank("REC::Traj")) recTraj = event.getBank("REC::Traj"); + + Particle recEl = null; + Particle recPr = null; + Particle recPip = null; + Particle recPim = null; + Particle recPro = null; + LorentzVector virtualPhoton = null; + LorentzVector hadronSystem = null; + LorentzVector virtualPhotonP = null; + LorentzVector hadronSystemP = null; + Particle beam = Particle.createWithPid(11, 0,0,ebeam, 0,0,0); + Particle target = Particle.createWithPid(2212, 0,0,0, 0,0,0); + if(event.hasBank("REC::Particle")==true && event.hasBank("REC::Track")){ + DataBank bank = event.getBank("REC::Particle"); + DataBank track = event.getBank("REC::Track"); + int rows = bank.rows(); + for(int loop = 0; loop < rows; loop++){ + int status = (int) Math.floor(Math.abs(bank.getShort("status", loop))/1000); /// + if(loop==0 && bank.getInt("pid", loop)==11 && status==2) { + recEl = new Particle( + bank.getInt("pid", loop), + bank.getFloat("px", loop), + bank.getFloat("py", loop), + bank.getFloat("pz", loop), + bank.getFloat("vx", loop), + bank.getFloat("vy", loop), + bank.getFloat("vz", loop)); + for(int j=0; j0){ + histos.get("General").getH2F("hi_q2w").fill(hadronSystem.mass(),-virtualPhoton.mass2()); + histos.get("General").getH1F("hi_w").fill(hadronSystem.mass()); + histos.get("General").getH2F("hi_z_phi").fill(Math.toDegrees(recEl.phi()), recEl.vz()); + histos.get("General").getH2F("hi_w_phi").fill(Math.toDegrees(recEl.phi()), hadronSystem.mass()); + histos.get("General").getH2F("hi_w_theta").fill(Math.toDegrees(recEl.theta()), hadronSystem.mass()); + histos.get("General").getH2F("hi_el").fill(recEl.p(),Math.toDegrees(recEl.theta())); + histos.get("Elastic").getH1F("hi_w").fill(hadronSystem.mass()); + histos.get("Elastic").getH2F("hi_w_phi").fill(Math.toDegrees(recEl.phi()), hadronSystem.mass()); + histos.get("Elastic").getH2F("hi_w_theta").fill(Math.toDegrees(recEl.theta()), hadronSystem.mass()); + histos.get("Elastic").getH2F("hi_w_p").fill(recEl.p(),hadronSystem.mass()); + histos.get("W").getH1F("hi_w_" + secEl).fill(hadronSystem.mass()); + histos.get("Vertex").getH1F("h1_vz_el").fill(recEl.vz()); + histos.get("Vertex").getH2F("h2_theta_el").fill(recEl.vz(), Math.toDegrees(recEl.theta())); + histos.get("Vertex").getH2F("h2_phi_el").fill(recEl.vz()), Math.toDegrees(recEl.phi()); + histos.get("Vertex").getH2F("h2_p_el").fill(recEl.vz(), recEl.p()); + histos.get("Vertex-Sector").getH1F("h1_vz_el_"+secEl).fill(recEl.vz()); + } + if(hadronSystem.mass()<1.1) { + if(recPr != null) { + double phPr = Math.toDegrees(recPr.phi()); + if(phPr < phEl) phPr +=360; + histos.get("Proton").getH2F("hi_pr").fill(recPr.p(),Math.toDegrees(recPr.theta())); + // dg_proton.getH2F("hi_phi").fill(Math.toDegrees(recEl.phi()),Math.toDegrees(recPr.phi())); + // dg_proton.getH2F("hi_dphi").fill(Math.toDegrees(recPr.phi()),phPr-phEl-180); + histos.get("Phi").getH1F("hi_dphi_" + secEl).fill(phPr-phEl); + histos.get("Proton").getH1F("hi_dphi").fill(phPr-phEl); + histos.get("Proton").getH2F("hi_phi_dphi").fill(Math.toDegrees(recPr.phi()),phPr-phEl); + histos.get("Proton").getH2F("hi_phi_dz").fill(Math.toDegrees(recPr.phi()),recPr.vz()-recEl.vz()); + histos.get("Proton").getH1F("hi_dz").fill(recPr.vz()-recEl.vz()); + if(Math.abs(phPr-phEl-180)<10 && recPr.getProperty("NDF")>=0) { + histos.get("Proton").getH2F("hi_p_dp").fill(recPr.p(),-0.93832+recPr.e()+recEl.p()); + histos.get("Proton").getH1F("hi_dp").fill(-0.93832+recPr.e()+recEl.p()); + histos.get("Proton").getH1F("hi_dtheta").fill(Math.toDegrees(recPr.theta())-histos.get("Proton").getF1D("f1_pr").evaluate(recPr.p())); + histos.get("Beam").getH1F("hi_beam_" + secEl).fill((-0.93832+recPr.e()+recEl.p())); + if(lundfile!=null) { + PhysicsEvent ev = writeToLund(ebeam, recEl); + if(ev!=null) lundfile << ev.toLundString(); + } + } + } + } + if(recPip!=null) { + histos.get("Vertex").getH1F("h1_vz_pi+").fill(recPip.vz()); + histos.get("Vertex").getH2F("h2_theta_pi+").fill(recPip.vz(), Math.toDegrees(recPip.theta())); + histos.get("Vertex").getH2F("h2_phi_pi+").fill(recPip.vz(), Math.toDegrees(recPip.phi())); + histos.get("Vertex").getH2F("h2_p_pi+").fill(recPip.vz(), recPip.p()); + histos.get("Vertex-Sector").getH1F("h1_vz_pi+_"+((int) recPip.getProperty("sector"))).fill(recPip.vz()); + } + if(recPim!=null) { + histos.get("Vertex").getH1F("h1_vz_pi-").fill(recPim.vz()); + histos.get("Vertex").getH2F("h2_theta_pi-").fill(recPim.vz(), Math.toDegrees(recPim.theta())); + histos.get("Vertex").getH2F("h2_phi_pi-").fill(recPim.vz(), Math.toDegrees(recPim.phi())); + histos.get("Vertex").getH2F("h2_p_pi-").fill(recPim.vz(), recPim.p()); + histos.get("Vertex-Sector").getH1F("h1_vz_pi-_"+((int) recPim.getProperty("sector"))).fill(recPim.vz()); + } + } + if(recEl!=null && recPip!=null && recPim!=null) { + recPro = new Particle(); + recPro.copy(target); + recPro.combine(beam, +1); + recPro.combine(recEl, -1); + recPro.combine(recPip, -1); + recPro.combine(recPim, -1); + Particle rho = new Particle(); + rho.copy(recPip); + rho.combine(recPim, +1); + histos.get("2pi").getH1F("mxt1").fill(recPro.mass()); + histos.get("2pi").getH1F("mit1").fill(rho.mass()); + } + else if(recEl!=null && recPip!=null && recPro!=null && recPim==null) { + recPim = new Particle(); + recPim.copy(target); + recPim.combine(beam, +1); + recPim.combine(recEl, -1); + recPim.combine(recPip, -1); + recPim.combine(recPro, -1); + Particle rho = new Particle(); + rho.copy(recPip); + rho.combine(recPim, +1); + histos.get("2pi").getH1F("mxt2").fill(recPim.mass2()); + histos.get("2pi").getH1F("mit2").fill(rho.mass()); + } + else if(recEl!=null && recPim!=null && recPro!=null && recPip==null) { + recPip = new Particle(); + recPip.copy(target); + recPip.combine(beam, +1); + recPip.combine(recEl, -1); + recPip.combine(recPim, -1); + recPip.combine(recPro, -1); + Particle rho = new Particle(); + rho.copy(recPip); + rho.combine(recPim, +1); + histos.get("2pi").getH1F("mxt3").fill(recPip.mass2()); + histos.get("2pi").getH1F("mit3").fill(rho.mass()); + } + else if(recEl!=null && recPip!=null && recPim==null && recPro==null) { + Particle neutron = new Particle(); + neutron.copy(target); + neutron.combine(beam, +1); + neutron.combine(recEl, -1); + neutron.combine(recPip, -1); + Particle W = new Particle(); + W.copy(target); + W.combine(beam, +1); + W.combine(recEl, -1); + histos.get("1pi").getH2F("W1").fill(Math.toDegrees(recEl.phi()), neutron.mass()); + histos.get("1pi").getH1F("mxt1").fill(neutron.mass()); + } + else if(recEl!=null && recPro!=null && recPim==null && recPip==null) { + Particle pizero = new Particle(); + pizero.copy(target); + pizero.combine(beam, +1); + pizero.combine(recEl, -1); + pizero.combine(recPro, -1); + Particle W = new Particle(); + W.copy(target); + W.combine(beam, +1); + W.combine(recEl, -1); + histos.get("1pi").getH2F("W2").fill(Math.toDegrees(recEl.phi()), pizero.mass2()); + histos.get("1pi").getH1F("mxt2").fill(pizero.mass2()); + } + } + + } + + public void analyzeHistos() { + Logger.getLogger("org.freehep.math.minuit").setLevel(Level.WARNING); + + fitW(histos.get("General").getH1F("hi_w"), 0.8, 1.1, 0.05); + fitW(histos.get("Elastic").getH1F("hi_w"), 0.8, 1.1, 0.05); + fitW(histos.get("2pi").getH1F("mxt1"), 0.8, 1.1, 0.05); + fitW(histos.get("2pi").getH1F("mxt2"), -0.1, 0.2, 0.05); + fitW(histos.get("2pi").getH1F("mxt3"), -0.1, 0.2, 0.05); + fitW(histos.get("1pi").getH1F("mxt1"), 0.8, 1.1, 0.05); + fitW(histos.get("1pi").getH1F("mxt2"), -0.1, 0.2, 0.05); + for(int sector=1; sector <= 6; sector++) { + fitW(histos.get("W").getH1F("hi_w_" + sector), 0.8, 1.1, 0.05); + fitGauss(histos.get("Phi").getH1F("hi_dphi_" + sector)); + fitGauss(histos.get("Beam").getH1F("hi_beam_" + sector)); + } + fitGauss(histos.get("Proton").getH1F("hi_dp")); + fitGauss(histos.get("Proton").getH1F("hi_dtheta")); + fitGauss(histos.get("Proton").getH1F("hi_dphi")); + fitGauss(histos.get("Proton").getH1F("hi_dz")); +// fitW(histos.get("Vertex").getH1F("h1_vz_el"), scWindow+targetPos-2,scWindow+targetPos+2, 1); +// fitW(histos.get("Vertex").getH1F("h1_vz_pi+"), scWindow+targetPos-2,scWindow+targetPos+2, 1); +// fitW(histos.get("Vertex").getH1F("h1_vz_pi-"), scWindow+targetPos-2,scWindow+targetPos+2, 1); + } + + public EmbeddedCanvasTabbed drawHistos(String optStats) { + + EmbeddedCanvasTabbed canvas = null; + + String[] titles = (String[]) histos.keySet().toArray(); + for(int i=0; i dsList = histos.get(key).getData(i); + for(IDataSet ds : dsList){ + // System.out.println("\t --> " + ds.getName()); + dir.addDataSet(ds); + } + } + } + dir.writeFile(filename); + } + + public void readHistos(String filename) { + this.createHistos(); + TDirectory dir = new TDirectory(); + dir.readFile(filename); + System.out.println(dir.getDirectoryList()); + for(String key : histos.keySet()) { + String folder = key + "/"; + System.out.println("Reading from: " + folder); + DataGroup group = this.histos.get(key); + int nrows = group.getRows(); + int ncols = group.getColumns(); + int nds = nrows*ncols; + DataGroup newGroup = new DataGroup(ncols,nrows); + for(int i = 0; i < nds; i++){ + List dsList = group.getData(i); + for(IDataSet ds : dsList){ + System.out.println("\t --> " + ds.getName()); + newGroup.addDataSet(dir.getObject(folder, ds.getName()),i); + } + } + this.histos.replace(key, newGroup); + } + } + + private void fitW(H1F hiw, double min, double max, double sigma) { + + F1D f1w = new F1D("f1_w", "[amp]*gaus(x,[mean],[sigma])", min, max); + // get histogram maximum in the rane 0.7-1.1 + int i1=hiw.getXaxis().getBin(min); + int i2=hiw.getXaxis().getBin(max); + double hiMax=0; + int imax=i1; + for(int i=i1; i<=i2; i++) { + if(hiMax inputList = null; + + Kinematics analysis = null; + + if(parser.getCommand().equals("-process")) { + int maxEvents = parser.getOptionParser("-process").getOption("-nevent").intValue(); + double beamEnergy = parser.getOptionParser("-process").getOption("-beam").doubleValue(); + String vertexPar = parser.getOptionParser("-process").getOption("-vertpar").stringValue(); + String namePrefix = parser.getOptionParser("-process").getOption("-o").stringValue(); + String histoName = "histo.hipo"; + if(!namePrefix.isEmpty()) { + histoName = namePrefix + "_" + histoName; + } + optStats = parser.getOptionParser("-process").getOption("-stats").stringValue(); + openWindow = parser.getOptionParser("-process").getOption("-display").intValue()!=0; + if(!openWindow) System.setProperty("java.awt.headless", "true"); + + inputList = parser.getOptionParser("-process").getInputList(); + if(inputList.isEmpty()==true){ + parser.printUsage(); + System.out.println("\n >>>> error: no input file is specified....\n"); + System.exit(0); + } + + analysis = new Kinematics(beamEnergy, vertexPar); + + ProgressPrintout progress = new ProgressPrintout(); + + int counter = -1; + analysis.createHistos(); + for(String inputFile : inputList){ + HipoDataSource reader = new HipoDataSource(); + reader.open(inputFile); + + + while (reader.hasEvent()) { + + counter++; + + DataEvent event = reader.getNextEvent(); + analysis.processEvent(event); + + progress.updateStatus(); + if(maxEvents>0){ + if(counter>=maxEvents) break; + } + } + progress.showStatus(); + reader.close(); + } + analysis.analyzeHistos(); + analysis.saveHistos(histoName); + } + + if(parser.getCommand().equals("-plot")) { + double beamEnergy = parser.getOptionParser("-plot").getOption("-beam").doubleValue(); + String vertexPar = parser.getOptionParser("-plot").getOption("-vertpar").stringValue(); + optStats = parser.getOptionParser("-plot").getOption("-stats").stringValue(); + openWindow = parser.getOptionParser("-plot").getOption("-display").intValue()!=0; + if(!openWindow) System.setProperty("java.awt.headless", "true"); + + inputList = parser.getOptionParser("-plot").getInputList(); + if(inputList.isEmpty()==true){ + parser.printUsage(); + System.out.println("\n >>>> error: no input file is specified....\n"); + System.exit(0); + } + + analysis = new Kinematics(beamEnergy, vertexPar); + analysis.readHistos(inputList.get(0)); + analysis.analyzeHistos(); + } + + + if(openWindow) { + JFrame frame = new JFrame(inputList.get(0)); + frame.setSize(1200, 800); + frame.add(analysis.drawHistos(optStats)); + frame.setLocationRelativeTo(null); + frame.setVisible(true); + } + } +} From 653e2c4b15f718a4c81c747e471ecb05d4185235 Mon Sep 17 00:00:00 2001 From: Raffaella De Vita Date: Tue, 27 Jan 2026 16:57:09 -0500 Subject: [PATCH 05/15] more logging on fit failures and improvements to the RG_L vertex fit --- .../java/org/clas/dc/alignment/Alignment.java | 21 ++++-- .../java/org/clas/dc/alignment/Constants.java | 1 + .../java/org/clas/dc/alignment/Histo.java | 67 +++++++++++-------- 3 files changed, 56 insertions(+), 33 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java index f387e8e..f8cc47b 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java @@ -207,6 +207,8 @@ else if(Math.abs(vertex)==9) Constants.initTargetPars(Constants.RGKSPRING2024); else if(Math.abs(vertex)==10) Constants.initTargetPars(Constants.RGESPRING2024); + else if(Math.abs(vertex)==11) + Constants.initTargetPars(Constants.RGLSPRING2025); } if(vertex<0) vertex=0; } @@ -256,18 +258,25 @@ private void printExclusionStats() { LOGGER.log(Level.WARNING, "\nEcluded measurements because of failed fits:"); for (int is = 0; is < Constants.NSECTOR; is++) { String si = ""; + String ki = ""; int nexclude = 0; for (int it = 1; it < thetaBins.length; it++) { for (int ip = 1; ip < phiBins.length; ip++) { for (int il = 0; il < Constants.NLAYER + Constants.NTARGET; il++) { if (Constants.MEASWEIGHTS[is][it][ip][il] == 0) { nexclude++; - si += String.format("\n\t\ttheta bin=%d phi bin=%d layer=%d", it, ip, il); + si += String.format("\n\t\ttheta bin=%d phi bin=%d layer=%d\t", it, ip, il); + for(String key : histos.keySet()) { + si += " " + (histos.get(key).getParStatus(is+1, it, ip)[il] ? 1 : 0 ); + } } } } } - si = "\tSector " + (is+1) + String.format(": %d/%d", nexclude, (thetaBins.length-1)*(phiBins.length-1)*(Constants.NLAYER+Constants.NTARGET)) + si; + for(String key : histos.keySet()) { + ki += String.format("%6s", key); + } + si = "\tSector " + (is+1) + String.format(": %d/%d", nexclude, (thetaBins.length-1)*(phiBins.length-1)*(Constants.NLAYER+Constants.NTARGET)) + "\t\t\t" + ki + si; LOGGER.log(Level.WARNING, si); } } @@ -805,7 +814,7 @@ private DataGroup getVertexGraph(String parameter, Table alignment) { gr_fit.addPoint(phi, shiftRes/Constants.SCALE, 0.0, errorRes/Constants.SCALE); } } - gr_fit.setTitle("Layer " + (il+1)); + gr_fit.setTitle("Layer " + il); gr_fit.setTitleX("#phi (deg)"); gr_fit.setTitleY("#Deltaz (cm)"); gr_fit.setMarkerColor(this.markerColor[it-1]); @@ -1363,9 +1372,9 @@ public static void main(String[] args){ if(parser.getCommand().equals("-analyze")) { String namePrefix = parser.getOptionParser("-analyze").getOption("-o").stringValue(); - String histoName = "histo.hipo"; + String histoName = null; if(!namePrefix.isEmpty()) { - histoName = namePrefix + "_" + histoName; + histoName = namePrefix + "_histo.hipo"; } String optStats = parser.getOptionParser("-analyze").getOption("-stats").stringValue(); int residuals = parser.getOptionParser("-analyze").getOption("-residuals").intValue(); @@ -1393,7 +1402,7 @@ public static void main(String[] args){ align.initConstants(11, initVar, previousVar, compareVar); align.readHistos(inputHisto, optStats); align.analyzeHistos(residuals, vertexFit, vertexPar, testFit); - align.saveHistos(histoName); + if(histoName != null) align.saveHistos(histoName); } if(parser.getCommand().equals("-fit")) { diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java index 52b9668..98ebdd4 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java @@ -98,6 +98,7 @@ public class Constants { public static final double[] RGDFALL2023 = {-3.0, 5.0, 3.0, 27.1}; public static final double[] RGKSPRING2024 = {-0.5, 5.0, 0.0, 26.5}; // warm public static final double[] RGESPRING2024 = {-5.5, 2.2, 2.0, 25.9}; + public static final double[] RGLSPRING2025 = {19.2,50.4, 0.9, 0.0}; // target parameters used for vertex fit initialization public static double TARGETPOS = DEFAULT[0]; diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 536b72c..1c42ad6 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -49,14 +49,15 @@ public class Histo { private DataGroup[][][] leftright = null; // indices are theta bin, phi bin and sector, datagroup is 6x6 and contains layers private DataGroup[][] vertex = null; // indices are theta bin, phi bin and sector, datagroup is 6x6 and contains sectors - private double[][][][] zeroes = null; - private double[][][][] timeValues = null; - private double[][][][] timeSigmas = null; - private double[][][][] lrValues = null; - private double[][][][] parValues = null; - private double[][][][] parErrors = null; - private double[][][][] parSigmas = null; - private double[][] beamOffset= {{0, 0}, {0, 0}}; + private double[][][][] zeroes = null; + private double[][][][] timeValues = null; + private double[][][][] timeSigmas = null; + private double[][][][] lrValues = null; + private double[][][][] parValues = null; + private double[][][][] parErrors = null; + private double[][][][] parSigmas = null; + private boolean[][][][] parStatus = null; + private double[][] beamOffset= {{0, 0}, {0, 0}}; private Bin[] thetaBins = null; private Bin[] phiBins = null; @@ -137,6 +138,7 @@ private void createHistos(String optStats) { this.parValues = new double[nSector][thetaBins.length][phiBins.length][nLayer+nTarget]; this.parErrors = new double[nSector][thetaBins.length][phiBins.length][nLayer+nTarget]; this.parSigmas = new double[nSector][thetaBins.length][phiBins.length][nLayer+nTarget]; + this.parStatus = new boolean[nSector][thetaBins.length][phiBins.length][nLayer+nTarget]; this.timeValues = new double[nSector][thetaBins.length][phiBins.length][nLayer+nTarget]; this.timeSigmas = new double[nSector][thetaBins.length][phiBins.length][nLayer+nTarget]; this.lrValues = new double[nSector][thetaBins.length][phiBins.length][nLayer+nTarget]; @@ -630,6 +632,7 @@ public void analyzeHisto(int fit, int vertexFit) { this.parValues[is][it][ip][il] = 0; this.parErrors[is][it][ip][il] = 0; this.parSigmas[is][it][ip][il] = 1; + this.parStatus[is][it][ip][il] = false; } for(int l=1; l<=nLayer; l++) { if(tres) { @@ -652,6 +655,7 @@ public void analyzeHisto(int fit, int vertexFit) { this.parValues[is][it][ip][l] = hres.getFunction().getParameter(1); this.parErrors[is][it][ip][l] = hres.getFunction().parameter(1).error(); this.parSigmas[is][it][ip][l] = hres.getFunction().getParameter(2); + this.parStatus[is][it][ip][l] = true; if(!shift) { this.parErrors[is][it][ip][l] = Math.max(this.parErrors[is][it][ip][l],(Constants.RESMAX-Constants.RESMIN)/Constants.RESBINS/2); // if(l>24) this.parErrors[is][it][ip][l] *= 2; @@ -668,6 +672,7 @@ public void analyzeHisto(int fit, int vertexFit) { this.parValues[is][it][ip][0] = hvtx.getFunction().getParameter(1)*Constants.SCALE; this.parErrors[is][it][ip][0] = hvtx.getFunction().parameter(1).error()*Constants.SCALE; this.parSigmas[is][it][ip][0] = hvtx.getFunction().getParameter(2)*Constants.SCALE; + this.parStatus[is][it][ip][0] = true; if(!shift) { this.parValues[is][it][ip][0] -= Constants.TARGETPOS*Constants.SCALE; this.parErrors[is][it][ip][0] = Math.max(this.parErrors[is][it][ip][0], Constants.SCALE*dx/2); @@ -683,6 +688,7 @@ public void analyzeHisto(int fit, int vertexFit) { this.parValues[is][it][ip][nLayer+nTarget-1] = (hvtx.getFunction().getParameter(iscw)-Constants.SCEXIT)*Constants.SCALE; this.parErrors[is][it][ip][nLayer+nTarget-1] = Math.max(hvtx.getFunction().parameter(iscw).error()*Constants.SCALE, Constants.SCALE*dx); this.parSigmas[is][it][ip][nLayer+nTarget-1] = hvtx.getFunction().getParameter(2)*Constants.SCALE; + this.parStatus[is][it][ip][nLayer+nTarget-1] = true; } else { Constants.MEASWEIGHTS[is][it][ip][nLayer+nTarget-1]=0; @@ -691,6 +697,7 @@ public void analyzeHisto(int fit, int vertexFit) { this.parValues[is][it][ip][nLayer+nTarget-2] = (hvtx.getFunction().getParameter(itl)-Constants.TARGETLENGTH)*Constants.SCALE; this.parErrors[is][it][ip][nLayer+nTarget-2] = Math.max(hvtx.getFunction().parameter(itl).error()*Constants.SCALE, Constants.SCALE*dx); this.parSigmas[is][it][ip][nLayer+nTarget-2] = hvtx.getFunction().getParameter(2)*Constants.SCALE; + this.parStatus[is][it][ip][nLayer+nTarget-2] = true; } else { Constants.MEASWEIGHTS[is][it][ip][nLayer+nTarget-2]=0; @@ -755,6 +762,16 @@ public double[][] getBeamOffset() { return this.beamOffset; } + public boolean[] getParStatus(int sector, int itheta, int iphi) { + if(sector<1 || sector>6) + throw new IllegalArgumentException("Error: invalid sector="+sector); + if(itheta<0 || itheta>=thetaBins.length) + throw new IllegalArgumentException("Error: invalid theta bin="+itheta); + if(iphi<0 || iphi>phiBins.length) + throw new IllegalArgumentException("Error: invalid phi bin="+iphi); + return this.parStatus[sector-1][itheta][iphi]; + } + public double[] getParValues(int sector, int itheta, int iphi) { return this.getParValues("", sector, itheta, iphi); } @@ -1447,7 +1464,6 @@ public static void fitRGDVertex(H1F histo) { } /** -<<<<<<< HEAD * 3-peaks vertex fitting function * Peaks correspond to: target windows and scattering chamber exit window * Initialized according to: @@ -1637,26 +1653,20 @@ public static void fitRGEVertex(H1F histo) { /** * 2-peaks vertex fitting function - * Peaks correspond to: target windows and scattering chamber exit window + * Peaks correspond to: target windows and beam pipe exit window * Initialized according to: * - chosen target length (TARGETLENGTH), * - target exit window position (TARGETPOS) - * Includes a wide Gaussian to account for target gas + * - distance between the beam pipe exit window an target entrance window (WINDOWDIST) + * Includes a wide Gaussian and a polynomial to account for target gas * @param histo */ public static void fitRGLVertex(H1F histo) { int nbin = histo.getData().length; double dx = histo.getDataX(1)-histo.getDataX(0); - //find downstream window + //find windows int ibin1 = Histo.getMaximumBinBetween(histo, histo.getDataX(0), Constants.TARGETCENTER); int ibin2 = Histo.getMaximumBinBetween(histo, Constants.TARGETCENTER, histo.getDataX(nbin-1)); - //check if the found maximum is the first or second peak, ibin is tentative upstream window -// int ibin1 = Math.max(0, ibin0 - (int)(Constants.TARGETLENGTH/dx)); - // int ibin2 = Math.min(nbin-1, ibin0 + (int)(Constants.TARGETLENGTH/dx)); - // if(histo.getBinContent(ibin1) Date: Thu, 29 Jan 2026 19:14:57 -0500 Subject: [PATCH 06/15] improved RG-L vertex fit, added histogram rebinning, better fit failures logging and fit test --- .../java/org/clas/dc/alignment/Alignment.java | 35 +++++----- .../java/org/clas/dc/alignment/Constants.java | 22 ++++--- .../java/org/clas/dc/alignment/Histo.java | 64 +++++++++++++++---- 3 files changed, 84 insertions(+), 37 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java index f8cc47b..d481716 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java @@ -239,13 +239,13 @@ public void addHistoSet(String name, Histo histo) { this.histos.put(name, histo); } - public void analyzeHistos(int resFit, int vertexFit, String vertexPar, boolean test) { + public void analyzeHistos(int resFit, int vertexFit, String vertexPar, String test) { this.printConfig(resFit, "resFit", ""); this.printConfig(vertexFit, "vertexFit", ""); this.initMeasurementWeights(1.0); this.initVertexPar(vertexFit, vertexPar); for(String key : histos.keySet()) { - if(test && !key.equals("nominal")) continue; + if(histos.containsKey(test) && !key.equals(test)) continue; LOGGER.log(LEVEL,"\nAnalyzing histos for variation " + key); histos.get(key).analyzeHisto(resFit, vertexFit); for(int i=0; i0 && pars.length<=DEFAULT.length) { TARGETPOS = pars[0]; TARGETCENTER = pars[0]; + System.out.println("[CONFIG] target parameters set to:"); + System.out.println(" - TARGETPOS = " + TARGETPOS); if(pars.length>1) { TARGETLENGTH = pars[1]; TARGETCENTER = pars[0]-pars[1]/2; + System.out.println(" - TARGETLENGTH = " + TARGETLENGTH); } - if(pars.length>2) + if(pars.length>2) { WINDOWDIST = pars[2]; - if(pars.length>3) + System.out.println(" - WINDOWDIST = " + WINDOWDIST); + } + if(pars.length>3) { SCEXIT = pars[3]; + System.out.println(" - SCEXIT = " + SCEXIT); + } } else { System.out.println("[WARNING] wrong number of target parameters. Number is " + pars.length + " instead of [1:4]"); } - System.out.println("[CONFIG] target parameters set to:"); - System.out.println(" - TARGETPOS = " + TARGETPOS); - System.out.println(" - TARGETLENGTH = " + TARGETLENGTH); - System.out.println(" - WINDOWDIST = " + WINDOWDIST); - System.out.println(" - SCEXIT = " + SCEXIT); } } diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 1c42ad6..1d48c48 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -668,6 +668,7 @@ public void analyzeHisto(int fit, int vertexFit) { } H1F hvtx = vertex[it][ip].getH1F("hi-S"+s); double dx = hvtx.getDataX(1)-hvtx.getDataX(0); + System.out.print(String.format("\tsector=%1d theta bin=%1d phi bin=%1d layer=%2d",s,it,ip,0)); if(Histo.fitVertex(vertexFit, hvtx)) { this.parValues[is][it][ip][0] = hvtx.getFunction().getParameter(1)*Constants.SCALE; this.parErrors[is][it][ip][0] = hvtx.getFunction().parameter(1).error()*Constants.SCALE; @@ -709,6 +710,7 @@ public void analyzeHisto(int fit, int vertexFit) { Constants.MEASWEIGHTS[is][it][ip][nLayer+nTarget-2]=0; Constants.MEASWEIGHTS[is][it][ip][nLayer+nTarget-1]=0; } + System.out.print("\r"); } } } @@ -1662,28 +1664,33 @@ public static void fitRGEVertex(H1F histo) { * @param histo */ public static void fitRGLVertex(H1F histo) { + + double threshold = 100; + if(histo.getMax() Date: Tue, 10 Feb 2026 12:54:25 -0500 Subject: [PATCH 07/15] more vertex fit improvements --- .../java/org/clas/dc/alignment/Alignment.java | 3 +- .../java/org/clas/dc/alignment/Constants.java | 6 +-- .../java/org/clas/dc/alignment/Histo.java | 54 +++++++++---------- 3 files changed, 32 insertions(+), 31 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java index d481716..9c35c16 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java @@ -822,7 +822,8 @@ private DataGroup getVertexGraph(String parameter, Table alignment) { gr_fit.setTitle("Layer " + il); gr_fit.setTitleX("#phi (deg)"); gr_fit.setTitleY("#Deltaz (cm)"); - gr_fit.setMarkerColor(this.markerColor[it-1]); + gr_fit.setMarkerColor(this.markerColor[(it-1)%this.markerColor.length]); + gr_fit.setMarkerStyle(this.markerStyle[2 + (int) ((it-1)/this.markerColor.length)]); gr_fit.setMarkerSize(this.markerSize); if(gr_fit.getDataSize(0)>0) { if(!graphs.containsKey(i)) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java index ddede3a..a4d547a 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java @@ -61,9 +61,9 @@ public class Constants { "r3_x", "r3_y", "r3_z", "r3_cx", "r3_cy", "r3_cz"}; // parameter step size: set to 0 to fix the parameter - public static double[] PARSTEP = { 0.2, 0.0, 0.0, 0, 0.0, 0.0, - 0.2, 0.0, 0.2, 0, 0.0, 0.0, - 0.2, 0.0, 0.2, 0, 0.0, 0.0}; + public static double[] PARSTEP = { 0.2, 0.2, 0.2, 0, 0.2, 0.0, + 0.2, 0.2, 0.2, 0, 0.2, 0.0, + 0.2, 0.2, 0.2, 0, 0.2, 0.0}; // parameter max value public static double[] PARMAX = { 1.5, 1.5, 1.5, 0.5, 0.5, 0.5, diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 1d48c48..92f374f 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -1679,41 +1679,41 @@ public static void fitRGLVertex(H1F histo) { double ampU = histo.getBinContent(ibin1); double ampD = histo.getBinContent(ibin2); double sigma = 0.5; - double bg = histo.getBinContent((ibin1+ibin2)/2); - String function = "[ampU]*gaus(x,[exw]-[tll],[sigmaU])+" - + "[ampU]*gaus(x,[exw]-[tll]-[wd],[sigmaU])+" - + "[ampD]*gaus(x,[exw],[sigmaD])+" - + "[bg]*landau(x,[bgmean],[bgsigma])+" -// + "[p0]*landau(x,[p1],[p2])"; + double bg = ibin1>0 ? histo.getBinContent((ibin1+ibin2)/2) : histo.getBinContent(0); + String function = "[ampD]*gaus(x,[exw],[sigmaD])+" + + "[ampD]*0.1*landau(x,[exw],[sigmaD])+" + "[p0]+[p1]*x+[p2]*x*x"; F1D f1_vtx = new F1D("f"+histo.getName(), function, -10, 10); f1_vtx.setLineColor(2); f1_vtx.setLineWidth(2); - f1_vtx.setOptStat("11111111111111111"); - f1_vtx.setParameter(0, ampU/2); + f1_vtx.setOptStat("1111111111"); + f1_vtx.setParameter(0, ampD); f1_vtx.setParameter(1, meanD); - f1_vtx.setParameter(2, meanD-meanU);//Constants.TARGETLENGTH); - f1_vtx.setParLimits(2, Constants.TARGETLENGTH*0.9, Constants.TARGETLENGTH*1.1); - f1_vtx.setParameter(3, sigma*2); - f1_vtx.setParameter(4, Constants.WINDOWDIST); - f1_vtx.setParLimits(4, Constants.WINDOWDIST*0.99, Constants.WINDOWDIST*1.01); - f1_vtx.setParameter(5, ampD); - f1_vtx.setParameter(6, sigma); - f1_vtx.setParameter(7, bg/2); - f1_vtx.setParLimits(7, 0, bg*2); - f1_vtx.setParameter(8, meanD);//-Constants.TARGETLENGTH*0.5); - f1_vtx.setParLimits(8, meanD-2*sigma,meanD+6*sigma); - f1_vtx.setParameter(9, sigma*2); - f1_vtx.setParLimits(9, 0, sigma*8); - f1_vtx.setParameter(10, bg); - f1_vtx.setParLimits(10, 0, bg*2); -// f1_vtx.setParameter(11, meanU); -// f1_vtx.setParameter(12, Constants.TARGETLENGTH/3); + f1_vtx.setParameter(2, sigma); +// f1_vtx.setParameter(3, ampD/20); +// f1_vtx.setParLimits(3, ampD/40,ampD/10); +// f1_vtx.setParameter(4, meanD); +// f1_vtx.setParLimits(4, meanD-2*sigma,meanD+6*sigma); +// f1_vtx.setParameter(5, sigma*2); +// f1_vtx.setParLimits(5, sigma, sigma*4); + if(ibin1>0) { + int np = f1_vtx.getParameterEstimate().length; + f1_vtx = new F1D("f"+histo.getName(), function + + "+[ampU]*gaus(x,[exw]-[tl],[sigmaU])" + + "+[ampU]*gaus(x,[exw]-[tl]-[wd],[sigmaU])", -10, 10); + f1_vtx.setOptStat("11111111111111"); + f1_vtx.setParameter(np+0, ampU/2); + f1_vtx.setParameter(np+1, meanD-meanU);//Constants.TARGETLENGTH); + f1_vtx.setParLimits(np+1, Constants.TARGETLENGTH*0.9, Constants.TARGETLENGTH*1.1); + f1_vtx.setParameter(np+2, sigma*2); + f1_vtx.setParameter(np+3, Constants.WINDOWDIST); + f1_vtx.setParLimits(np+3, Constants.WINDOWDIST*0.99, Constants.WINDOWDIST*1.01); + } f1_vtx.setRange(Math.max(meanU-8*sigma,histo.getDataX(0)), Math.min(meanD+8*sigma,histo.getDataX(nbin-1))); histo.setFunction(f1_vtx); DataFitter.fit(f1_vtx, histo, "Q"); //No options uses error for sigma - if(!f1_vtx.isFitValid()) { + if(!f1_vtx.isFitValid() && ibin1>0) { meanU = f1_vtx.getParameter(1)-f1_vtx.getParameter(2); meanD = f1_vtx.getParameter(8); f1_vtx.setRange(Math.max(meanU-8*sigma,histo.getDataX(0)), @@ -1824,7 +1824,7 @@ public static int getMaximumBinBetween(H1F histo, double min, double max) { double x_val; double y_max_temp; double y_max = 0; - int max_bin_num = histo.getMaximumBin(); + int max_bin_num = 0; for (int i = 0; i < nbin; i++) { x_val_temp = histo.getAxis().getBinCenter(i); if (x_val_temp >= min && x_val_temp <= max) { From 4ea2a5478377d6d79b57263e3d55be0d9ce7bdf2 Mon Sep 17 00:00:00 2001 From: Raffaella De Vita Date: Wed, 11 Feb 2026 18:11:07 -0500 Subject: [PATCH 08/15] more general plots on residuals quality and recalculation of fit residuals --- .../java/org/clas/dc/alignment/Histo.java | 101 ++++++++++++------ 1 file changed, 69 insertions(+), 32 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 92f374f..e4bf9ad 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -42,8 +42,9 @@ public class Histo { private DataGroup offset = null; private DataGroup calib = null; private DataGroup alpha = null; + private DataGroup track = null; private DataGroup doca = null; - private DataGroup[] wires = null; + private DataGroup[][] wires = null; private DataGroup[][][] residuals = null; // indices are theta bin, phi bin and sector, datagroup is 6x6 and contains layers private DataGroup[][][] time = null; // indices are theta bin, phi bin and sector, datagroup is 6x6 and contains layers private DataGroup[][][] leftright = null; // indices are theta bin, phi bin and sector, datagroup is 6x6 and contains layers @@ -128,8 +129,8 @@ private void createHistos(String optStats) { + (phiBins.length-1) + " phi bins " + "for variation " + name); this.residuals = new DataGroup[nSector][thetaBins.length][phiBins.length]; + this.wires = new DataGroup[nSector][2]; if(tres) { - this.wires = new DataGroup[nSector]; this.time = new DataGroup[nSector][thetaBins.length][phiBins.length]; this.leftright = new DataGroup[nSector][thetaBins.length][phiBins.length]; } @@ -158,6 +159,7 @@ private void createHistos(String optStats) { this.calib = new DataGroup(nSector, nSLayer); this.alpha = new DataGroup(nSector, nSLayer); + this.track = new DataGroup(nSector, nSLayer); this.doca = new DataGroup(nSector, nSLayer); for(int is=0; is getHits(Event event, Electron e) { int layer = hitBank.getInt("layer", i) + 6 * (superlayer - 1); int wire = hitBank.getInt("wire", i); int status = hitBank.getInt("status", i); + residual = doca*1E4/Math.cos(Math.toRadians(alpha))-(doca*lr*1E4-residual); Hit hit = new Hit(sector, layer, wire, residual, time, doca, alpha, lr, status); if(!allhits.containsKey(superlayer)) allhits.put(superlayer, new ArrayList<>()); @@ -862,21 +874,30 @@ public EmbeddedCanvasTabbed getElectronPlots() { } public EmbeddedCanvasTabbed plotHistos() { - EmbeddedCanvasTabbed canvas = new EmbeddedCanvasTabbed("Calibration", "Alpha", "Doca"); + EmbeddedCanvasTabbed canvas = new EmbeddedCanvasTabbed("Calibration", "Alpha", "Track", "Doca"); canvas.getCanvas("Calibration").draw(calib); canvas.getCanvas("Alpha").draw(alpha); for(EmbeddedPad pad : canvas.getCanvas("Alpha").getCanvasPads()) pad.getAxisZ().setLog(true); + canvas.getCanvas("Track").draw(track); canvas.getCanvas("Doca").draw(doca); + for(int is=0; is Date: Thu, 12 Feb 2026 13:57:58 -0500 Subject: [PATCH 09/15] bug fixes --- .../java/org/clas/dc/alignment/Alignment.java | 2 +- .../java/org/clas/dc/alignment/Histo.java | 43 +++++++++++-------- 2 files changed, 25 insertions(+), 20 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java index 9c35c16..eeb7e5a 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Alignment.java @@ -1239,7 +1239,7 @@ public static void main(String[] args){ parser.getOptionParser("-process").addOption("-frame" , "0", "translations defined in the CLAS12 tilted sector frame (0) or sector frame (1)"); parser.getOptionParser("-process").addOption("-global" , "0", "r1 translations defined as relative (0) or global (1) translations"); parser.getOptionParser("-process").addOption("-verbose" , "0", "global fit verbosity (1/0 = on/off)"); - parser.getOptionParser("-process").addOption("-test" , "nominal", "geometry variation to be analyze only for fit testing"); + parser.getOptionParser("-process").addOption("-test" , "", "geometry variation to be analyze only for fit testing"); // valid options for histogram-base analysis parser.addCommand("-analyze", "analyze histogram files"); diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index e4bf9ad..5cbeaf2 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -129,7 +129,7 @@ private void createHistos(String optStats) { + (phiBins.length-1) + " phi bins " + "for variation " + name); this.residuals = new DataGroup[nSector][thetaBins.length][phiBins.length]; - this.wires = new DataGroup[nSector][2]; + this.wires = new DataGroup[2][nSector]; if(tres) { this.time = new DataGroup[nSector][thetaBins.length][phiBins.length]; this.leftright = new DataGroup[nSector][thetaBins.length][phiBins.length]; @@ -201,13 +201,13 @@ private void createHistos(String optStats) { for(int ires=0; ires<2; ires++) { for(int is=0; is Date: Thu, 12 Feb 2026 15:03:17 -0500 Subject: [PATCH 10/15] histo reorganization --- .../java/org/clas/dc/alignment/Histo.java | 339 +++++++++--------- 1 file changed, 171 insertions(+), 168 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 5cbeaf2..fd3296c 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -70,7 +70,7 @@ public class Histo { private List shiftedFiles = null; private SchemaFactory schema = null; - private boolean tres = false; + private boolean extra = false; private boolean shift = false; private static final Logger LOGGER = Logger.getLogger(Constants.LOGGERNAME); @@ -85,13 +85,13 @@ public Histo(String name, List files, Bin[] thetabins, Bin[] phibins, do this.createHistos(optstats); } - public Histo(String name, List files, Bin[] thetabins, Bin[] phibins, boolean time, double[] vertexrange, String optstats) { + public Histo(String name, List files, Bin[] thetabins, Bin[] phibins, boolean extraplots, double[] vertexrange, String optstats) { this.name = name; this.thetaBins = thetabins; this.phiBins = phibins; this.minVtx = vertexrange[0]; this.maxVtx = vertexrange[1]; - this.tres = time; + this.extra = extraplots; this.nominalFiles = files; this.createHistos(optstats); } @@ -114,12 +114,12 @@ public Histo(String name, boolean shift, Bin[] thetabins, Bin[] phibins, double[ this.createHistos(optstats); } - public Histo(String name, boolean shift, Bin[] thetabins, Bin[] phibins, boolean time, double[] vertexrange, String optstats) { + public Histo(String name, boolean shift, Bin[] thetabins, Bin[] phibins, boolean extraplots, double[] vertexrange, String optstats) { this.name = name; this.thetaBins = thetabins; this.phiBins = phibins; this.shift = shift; - this.tres = time; + this.extra = extraplots; this.createHistos(optstats); } @@ -129,8 +129,7 @@ private void createHistos(String optStats) { + (phiBins.length-1) + " phi bins " + "for variation " + name); this.residuals = new DataGroup[nSector][thetaBins.length][phiBins.length]; - this.wires = new DataGroup[2][nSector]; - if(tres) { + if(extra) { this.time = new DataGroup[nSector][thetaBins.length][phiBins.length]; this.leftright = new DataGroup[nSector][thetaBins.length][phiBins.length]; } @@ -157,16 +156,54 @@ private void createHistos(String optStats) { maxVtx = Constants.VDFMAX; } - this.calib = new DataGroup(nSector, nSLayer); - this.alpha = new DataGroup(nSector, nSLayer); - this.track = new DataGroup(nSector, nSLayer); - this.doca = new DataGroup(nSector, nSLayer); + if(extra) { + this.calib = new DataGroup(nSector, nSLayer); + this.alpha = new DataGroup(nSector, nSLayer); + this.track = new DataGroup(nSector, nSLayer); + this.doca = new DataGroup(nSector, nSLayer); + this.wires = new DataGroup[2][nSector]; + for(int ires=0; ires<2; ires++) { + for(int is=0; is0) this.leftright[sector-1][it][ip].getH1F("hi-rL" + hit.layer).fill(hit.time); @@ -567,7 +564,7 @@ private List getHits(Event event, Electron e) { int layer = hitBank.getInt("layer", i) + 6 * (superlayer - 1); int wire = hitBank.getInt("wire", i); int status = hitBank.getInt("status", i); - residual = doca*1E4/Math.cos(Math.toRadians(alpha))-(doca*lr*1E4-residual); + residual = doca*lr*1E4/Math.cos(Math.toRadians(alpha))-(doca*lr*1E4-residual); Hit hit = new Hit(sector, layer, wire, residual, time, doca, alpha, lr, status); if(!allhits.containsKey(superlayer)) allhits.put(superlayer, new ArrayList<>()); @@ -633,10 +630,12 @@ public void analyzeHisto(int fit, int vertexFit) { Histo.fitVertex(vertexFit,electron.getH1F("hi-vtx")); for(int is=0; is Date: Fri, 13 Feb 2026 16:25:03 -0500 Subject: [PATCH 11/15] removing residual recalculation --- .../dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java | 1 - 1 file changed, 1 deletion(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index fd3296c..cdbc219 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -564,7 +564,6 @@ private List getHits(Event event, Electron e) { int layer = hitBank.getInt("layer", i) + 6 * (superlayer - 1); int wire = hitBank.getInt("wire", i); int status = hitBank.getInt("status", i); - residual = doca*lr*1E4/Math.cos(Math.toRadians(alpha))-(doca*lr*1E4-residual); Hit hit = new Hit(sector, layer, wire, residual, time, doca, alpha, lr, status); if(!allhits.containsKey(superlayer)) allhits.put(superlayer, new ArrayList<>()); From bb9e239b5f4562739024a71a352fc160683788c1 Mon Sep 17 00:00:00 2001 From: Raffaella De Vita Date: Fri, 13 Feb 2026 17:46:58 -0500 Subject: [PATCH 12/15] bug fix --- .../src/main/java/org/clas/dc/alignment/Histo.java | 13 ++++++------- 1 file changed, 6 insertions(+), 7 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index cdbc219..5698880 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -1713,14 +1713,8 @@ public static void fitRGLVertex(H1F histo) { f1_vtx.setParameter(0, ampD); f1_vtx.setParameter(1, meanD); f1_vtx.setParameter(2, sigma); -// f1_vtx.setParameter(3, ampD/20); -// f1_vtx.setParLimits(3, ampD/40,ampD/10); -// f1_vtx.setParameter(4, meanD); -// f1_vtx.setParLimits(4, meanD-2*sigma,meanD+6*sigma); -// f1_vtx.setParameter(5, sigma*2); -// f1_vtx.setParLimits(5, sigma, sigma*4); if(ibin1>0) { - int np = f1_vtx.getParameterEstimate().length; + int np = f1_vtx.getNPars(); f1_vtx = new F1D("f"+histo.getName(), function + "+[ampU]*gaus(x,[exw]-[tl],[sigmaU])" + "+[ampU]*gaus(x,[exw]-[tl]-[wd],[sigmaU])", -10, 10); @@ -1732,6 +1726,11 @@ public static void fitRGLVertex(H1F histo) { f1_vtx.setParameter(np+3, Constants.WINDOWDIST); f1_vtx.setParLimits(np+3, Constants.WINDOWDIST*0.99, Constants.WINDOWDIST*1.01); } + f1_vtx.setLineColor(2); + f1_vtx.setLineWidth(2); + f1_vtx.setParameter(0, ampD); + f1_vtx.setParameter(1, meanD); + f1_vtx.setParameter(2, sigma); f1_vtx.setRange(Math.max(meanU-8*sigma,histo.getDataX(0)), Math.min(meanD+8*sigma,histo.getDataX(nbin-1))); histo.setFunction(f1_vtx); From e8fa5240b5e74748ab92ea36598a683fa0800ec7 Mon Sep 17 00:00:00 2001 From: Raffaella De Vita Date: Thu, 26 Feb 2026 17:23:24 -0500 Subject: [PATCH 13/15] more vertex fit improvements --- .../java/org/clas/dc/alignment/Histo.java | 25 +++++++++---------- 1 file changed, 12 insertions(+), 13 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 5698880..660dcd7 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -1704,20 +1704,16 @@ public static void fitRGLVertex(H1F histo) { double sigma = 0.5; double bg = ibin1>0 ? histo.getBinContent((ibin1+ibin2)/2) : histo.getBinContent(0); String function = "[ampD]*gaus(x,[exw],[sigmaD])+" - + "[ampD]*0.1*landau(x,[exw],[sigmaD])+" + + "[ampDL]*landau(x,[exw],[sigmaD]*3)+" + "[p0]+[p1]*x+[p2]*x*x"; F1D f1_vtx = new F1D("f"+histo.getName(), function, -10, 10); - f1_vtx.setLineColor(2); - f1_vtx.setLineWidth(2); f1_vtx.setOptStat("1111111111"); - f1_vtx.setParameter(0, ampD); - f1_vtx.setParameter(1, meanD); - f1_vtx.setParameter(2, sigma); + int np = f1_vtx.getNPars(); if(ibin1>0) { - int np = f1_vtx.getNPars(); f1_vtx = new F1D("f"+histo.getName(), function + "+[ampU]*gaus(x,[exw]-[tl],[sigmaU])" - + "+[ampU]*gaus(x,[exw]-[tl]-[wd],[sigmaU])", -10, 10); + + "+[ampU]*gaus(x,[exw]-[tl]-[wd],[sigmaU])" + + "+[ampUL]*landau(x,[exw]-[tl]-[wd]*0.5,[sigmaU]*3)", -10, 10); f1_vtx.setOptStat("11111111111111"); f1_vtx.setParameter(np+0, ampU/2); f1_vtx.setParameter(np+1, meanD-meanU);//Constants.TARGETLENGTH); @@ -1725,20 +1721,23 @@ public static void fitRGLVertex(H1F histo) { f1_vtx.setParameter(np+2, sigma*2); f1_vtx.setParameter(np+3, Constants.WINDOWDIST); f1_vtx.setParLimits(np+3, Constants.WINDOWDIST*0.99, Constants.WINDOWDIST*1.01); + f1_vtx.setParameter(np+4, ampU*0.1); } f1_vtx.setLineColor(2); f1_vtx.setLineWidth(2); f1_vtx.setParameter(0, ampD); f1_vtx.setParameter(1, meanD); f1_vtx.setParameter(2, sigma); - f1_vtx.setRange(Math.max(meanU-8*sigma,histo.getDataX(0)), + f1_vtx.setParameter(3, ampD*0.1); + f1_vtx.setRange(Math.max(meanU-12*sigma,histo.getDataX(0)), Math.min(meanD+8*sigma,histo.getDataX(nbin-1))); histo.setFunction(f1_vtx); DataFitter.fit(f1_vtx, histo, "Q"); //No options uses error for sigma - if(!f1_vtx.isFitValid() && ibin1>0) { - meanU = f1_vtx.getParameter(1)-f1_vtx.getParameter(2); - meanD = f1_vtx.getParameter(8); - f1_vtx.setRange(Math.max(meanU-8*sigma,histo.getDataX(0)), + if((!f1_vtx.isFitValid() ||f1_vtx.getChiSquare()/f1_vtx.getNDF()>Constants.CHI2MAX) && ibin1>0) { + meanU = f1_vtx.getParameter(1)-f1_vtx.getParameter(np+1); + meanD = f1_vtx.getParameter(1); + sigma = f1_vtx.getParameter(np+2); + f1_vtx.setRange(Math.max(meanU-6*sigma,histo.getDataX(0)), Math.min(meanD+8*sigma,histo.getDataX(nbin-1))); DataFitter.fit(f1_vtx, histo, "Q"); //No options uses error for sigma } From 3b78fc773d149555fec603ecbb7aa2a51a1df840 Mon Sep 17 00:00:00 2001 From: Raffaella De Vita Date: Mon, 2 Mar 2026 18:51:50 -0500 Subject: [PATCH 14/15] finalized rgl vertex fits --- .../java/org/clas/dc/alignment/Constants.java | 8 ++--- .../java/org/clas/dc/alignment/Histo.java | 32 +++++++++++-------- 2 files changed, 22 insertions(+), 18 deletions(-) diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java index a4d547a..ddfddf8 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Constants.java @@ -42,7 +42,7 @@ public class Constants { public static int VDFBINS = 200; public static double VDFMIN = -5.0; public static double VDFMAX = 5.0; - public static double CHI2MAX = 2.5; + public static double CHI2MAX = 5; // global fit @@ -61,9 +61,9 @@ public class Constants { "r3_x", "r3_y", "r3_z", "r3_cx", "r3_cy", "r3_cz"}; // parameter step size: set to 0 to fix the parameter - public static double[] PARSTEP = { 0.2, 0.2, 0.2, 0, 0.2, 0.0, - 0.2, 0.2, 0.2, 0, 0.2, 0.0, - 0.2, 0.2, 0.2, 0, 0.2, 0.0}; + public static double[] PARSTEP = { 0.2, 0.2, 0.2, 0, 0.2, 0.2, + 0.2, 0.2, 0.2, 0, 0.2, 0.2, + 0.2, 0.2, 0.2, 0, 0.2, 0.2}; // parameter max value public static double[] PARMAX = { 1.5, 1.5, 1.5, 0.5, 0.5, 0.5, diff --git a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java index 660dcd7..ccf4425 100644 --- a/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java +++ b/dc/java/dc-alignment/src/main/java/org/clas/dc/alignment/Histo.java @@ -1694,8 +1694,8 @@ public static void fitRGLVertex(H1F histo) { int nbin = histo.getData().length; //find windows - int ibin1 = Histo.getMaximumBinBetween(histo, histo.getDataX(0), Constants.TARGETCENTER-Constants.TARGETLENGTH/2); - int ibin2 = Histo.getMaximumBinBetween(histo, Constants.TARGETCENTER+Constants.TARGETLENGTH/2, histo.getDataX(nbin-1)); + int ibin1 = Histo.getMaximumBinBetween(histo, histo.getDataX(0), Constants.TARGETCENTER-Constants.TARGETLENGTH/4); + int ibin2 = Histo.getMaximumBinBetween(histo, Constants.TARGETCENTER+Constants.TARGETLENGTH/4, histo.getDataX(nbin-1)); double meanU = histo.getDataX(ibin1); double meanD = histo.getDataX(ibin2); @@ -1704,24 +1704,26 @@ public static void fitRGLVertex(H1F histo) { double sigma = 0.5; double bg = ibin1>0 ? histo.getBinContent((ibin1+ibin2)/2) : histo.getBinContent(0); String function = "[ampD]*gaus(x,[exw],[sigmaD])+" - + "[ampDL]*landau(x,[exw],[sigmaD]*3)+" - + "[p0]+[p1]*x+[p2]*x*x"; + + "[ampDL]*landau(x,[exw],[sigmaD]*2.5)+" + + "[f0]+([p0]+[p1]*x+[p2]*x*x)"; F1D f1_vtx = new F1D("f"+histo.getName(), function, -10, 10); f1_vtx.setOptStat("1111111111"); int np = f1_vtx.getNPars(); if(ibin1>0) { f1_vtx = new F1D("f"+histo.getName(), function + + "*erf(x,[exw]-[tl]-0.9,-[sigmaUL])" + "+[ampU]*gaus(x,[exw]-[tl],[sigmaU])" - + "+[ampU]*gaus(x,[exw]-[tl]-[wd],[sigmaU])" - + "+[ampUL]*landau(x,[exw]-[tl]-[wd]*0.5,[sigmaU]*3)", -10, 10); + + "+[ampU]*gaus(x,[exw]-[tl]-0.9,[sigmaU])", -10, 10); f1_vtx.setOptStat("11111111111111"); - f1_vtx.setParameter(np+0, ampU/2); - f1_vtx.setParameter(np+1, meanD-meanU);//Constants.TARGETLENGTH); - f1_vtx.setParLimits(np+1, Constants.TARGETLENGTH*0.9, Constants.TARGETLENGTH*1.1); - f1_vtx.setParameter(np+2, sigma*2); - f1_vtx.setParameter(np+3, Constants.WINDOWDIST); - f1_vtx.setParLimits(np+3, Constants.WINDOWDIST*0.99, Constants.WINDOWDIST*1.01); - f1_vtx.setParameter(np+4, ampU*0.1); + f1_vtx.setParameter(np+0, meanD-meanU-Constants.WINDOWDIST/2);//Constants.TARGETLENGTH); + f1_vtx.setParLimits(np+0, Constants.TARGETLENGTH*0.9, Constants.TARGETLENGTH*1.1); +// f1_vtx.setParameter(np+1, Constants.WINDOWDIST); +// f1_vtx.setParLimits(np+1, Constants.WINDOWDIST*0.99, Constants.WINDOWDIST*1.01); + f1_vtx.setParameter(np+1, sigma*5); + f1_vtx.setParStep(np+1, 0); + f1_vtx.setParameter(np+2, ampU/2); + f1_vtx.setParameter(np+3, sigma*2); +// f1_vtx.setParameter(np+4, ampU*0.1); } f1_vtx.setLineColor(2); f1_vtx.setLineWidth(2); @@ -1729,11 +1731,13 @@ public static void fitRGLVertex(H1F histo) { f1_vtx.setParameter(1, meanD); f1_vtx.setParameter(2, sigma); f1_vtx.setParameter(3, ampD*0.1); + f1_vtx.setParameter(4, 0); + f1_vtx.setParameter(5, bg/4); f1_vtx.setRange(Math.max(meanU-12*sigma,histo.getDataX(0)), Math.min(meanD+8*sigma,histo.getDataX(nbin-1))); histo.setFunction(f1_vtx); DataFitter.fit(f1_vtx, histo, "Q"); //No options uses error for sigma - if((!f1_vtx.isFitValid() ||f1_vtx.getChiSquare()/f1_vtx.getNDF()>Constants.CHI2MAX) && ibin1>0) { + if((false && !f1_vtx.isFitValid() ||f1_vtx.getChiSquare()/f1_vtx.getNDF()>Constants.CHI2MAX) && ibin1>0) { meanU = f1_vtx.getParameter(1)-f1_vtx.getParameter(np+1); meanD = f1_vtx.getParameter(1); sigma = f1_vtx.getParameter(np+2); From dd1cbe32d371b6ce60752500bf7c736202ffb286 Mon Sep 17 00:00:00 2001 From: Raffaella De Vita Date: Thu, 26 Mar 2026 18:58:38 -0400 Subject: [PATCH 15/15] version bump --- dc/java/dc-alignment/pom.xml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/dc/java/dc-alignment/pom.xml b/dc/java/dc-alignment/pom.xml index 3807c0a..6835c00 100644 --- a/dc/java/dc-alignment/pom.xml +++ b/dc/java/dc-alignment/pom.xml @@ -3,7 +3,7 @@ 4.0.0 org.clas.detector dc-alignment - 2.4 + 2.5 jar