diff --git a/fcl/detsim/detsim_2d_icarus_SConly.fcl b/fcl/detsim/detsim_2d_icarus_SConly.fcl new file mode 100644 index 000000000..83af1ead0 --- /dev/null +++ b/fcl/detsim/detsim_2d_icarus_SConly.fcl @@ -0,0 +1,9 @@ +#include "detsim_2d_icarus.fcl" + +physics.producers.tpcsim: @local::icarus_simwire_wirecell_SConly + +physics.simulate: ["rns", "tpcsim"] +physics.stream: [] +outputs: {} + +process_name: TPCSim diff --git a/fcl/detsim/detsim_2d_icarus_fitFR.fcl b/fcl/detsim/detsim_2d_icarus_fitFR.fcl index 9f69921f2..7e9ceb549 100644 --- a/fcl/detsim/detsim_2d_icarus_fitFR.fcl +++ b/fcl/detsim/detsim_2d_icarus_fitFR.fcl @@ -1,4 +1,4 @@ #include "detsim_2d_icarus.fcl" -physics.producers.daq: @local::icarus_simwire_wirecell_fitSR physics.producers.daq.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-multitpc-sim-drift-simchannel.jsonnet"] +physics.producers.daq: @local::icarus_simwire_wirecell_fitSR_P0nom diff --git a/fcl/reco/Definitions/sp_filter_parameters.fcl b/fcl/reco/Definitions/sp_filter_parameters.fcl new file mode 100644 index 000000000..5929eeb91 --- /dev/null +++ b/fcl/reco/Definitions/sp_filter_parameters.fcl @@ -0,0 +1,15 @@ +BEGIN_PROLOG +Wiener_tight_U_sigma: 0.065 +Wiener_tight_U_power: 4.4 +Wiener_tight_V_sigma: 0.065 +Wiener_tight_V_power: 2.6 +Wiener_tight_W_sigma: 0.07 +Wiener_tight_W_power: 3.4 + +Wire_ind_sigma: 0.4 +Wire_col_sigma: 2.2 + +Gaus_wide_sigma: 0.06 + +gain_ADC_per_e: 0.01214 # ADC / e-. For frame_scale +END_PROLOG diff --git a/fcl/reco/Definitions/stage0_icarus_defs.fcl b/fcl/reco/Definitions/stage0_icarus_defs.fcl index c22d2f9a9..2a6f08f2a 100644 --- a/fcl/reco/Definitions/stage0_icarus_defs.fcl +++ b/fcl/reco/Definitions/stage0_icarus_defs.fcl @@ -19,12 +19,17 @@ #include "wcls-decode-to-sig-base.fcl" #include "icarus_FilterDataIntegrity.fcl" #include "icarus_FilterCRTPMTMatching.fcl" +#include "spana_icarus.fcl" +#include "sp_filter_parameters.fcl" BEGIN_PROLOG ### Analyzers employed during stage 0 processing ### icarus_stage0_analyzers: { + spanaE: @local::spana1d_east + spanaW: @local::spana1d_west + purityinfoana0: { module_type: "TPCPurityInfoAna" PurityInfoLabel: "purityana0" PrintInfo: false @@ -205,15 +210,19 @@ icarus_stage0_multiTPC_TPC: [ decon1droi, roifinder1d ] -icarus_stage0_multiTPC_2d_TPC: [ decon1droi, - roifinder1d, +icarus_stage0_multiTPC_2d_TPC_E:[ decon2droiEE, - decon2droiEW, + decon2droiEW + ] + +icarus_stage0_multiTPC_2d_TPC_W:[ decon2droiWE, - decon2droiWW, - roifinder2d + decon2droiWW ] +icarus_stage0_multiTPC_2d_TPC: [@sequence::icarus_stage0_multiTPC_2d_TPC_E, @sequence::icarus_stage0_multiTPC_2d_TPC_W] + + icarus_stage0_EastHits_TPC: [ gaushit1dTPCEW, gaushit1dTPCEE ] @@ -273,13 +282,14 @@ icarus_stage0_multiTPC: [ @sequence::icarus_stage0_multiTPC_TPC, ] icarus_stage0_2d_multiTPC: [ @sequence::icarus_stage0_multiTPC_2d_TPC, - @sequence::icarus_stage0_EastHits_TPC, - @sequence::icarus_stage0_WestHits_TPC, @sequence::icarus_stage0_EastHits2d_TPC, @sequence::icarus_stage0_WestHits2d_TPC, @sequence::icarus_purity_monitor ] - + +icarus_stage0_2d_multiTPC_E: [ @sequence::icarus_stage0_multiTPC_2d_TPC_E, @sequence::icarus_stage0_EastHits2d_TPC] +icarus_stage0_2d_multiTPC_W: [ @sequence::icarus_stage0_multiTPC_2d_TPC_W, @sequence::icarus_stage0_WestHits2d_TPC] + icarus_stage0_CRT: [ daqCRT, crthit, @@ -346,30 +356,128 @@ icarus_stage0_producers.decon2droi.wcls_main.outputers: icarus_stage0_producers.decon2droi.wcls_main.params.raw_input_label: "daqTPC" icarus_stage0_producers.decon2droi.wcls_main.params.tpc_volume_label: 0 icarus_stage0_producers.decon2droi.wcls_main.params.signal_output_form: "dense" +icarus_stage0_producers.decon2droi.wcls_main.logsinks: ["stdout"] +icarus_stage0_producers.decon2droi.wcls_main.loglevels: ["debug", "pgraph:info"] +icarus_stage0_producers.decon2droi.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced.jsonnet"] icarus_stage0_producers.decon2droiEE.wcls_main.inputers: ["wclsRawFrameSource:rfsrc0"] icarus_stage0_producers.decon2droiEE.wcls_main.outputers: ["wclsFrameSaver:spsaver0"] icarus_stage0_producers.decon2droiEE.wcls_main.params.raw_input_label: "daqTPCROI:PHYSCRATEDATATPCEE" icarus_stage0_producers.decon2droiEE.wcls_main.params.tpc_volume_label: 0 icarus_stage0_producers.decon2droiEE.wcls_main.params.signal_output_form: "dense" +icarus_stage0_producers.decon2droiEE.wcls_main.logsinks: ["stdout"] +icarus_stage0_producers.decon2droiEE.wcls_main.loglevels: ["debug", "pgraph:info"] +icarus_stage0_producers.decon2droiEE.wcls_main.plugins: [@sequence::icarus_stage0_producers.decon2droiEE.wcls_main.plugins, "WireCellHio", "WireCellPytorch"] +icarus_stage0_producers.decon2droiEE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced.jsonnet"] icarus_stage0_producers.decon2droiEW.wcls_main.inputers: ["wclsRawFrameSource:rfsrc1"] icarus_stage0_producers.decon2droiEW.wcls_main.outputers: ["wclsFrameSaver:spsaver1"] icarus_stage0_producers.decon2droiEW.wcls_main.params.raw_input_label: "daqTPCROI:PHYSCRATEDATATPCEW" icarus_stage0_producers.decon2droiEW.wcls_main.params.tpc_volume_label: 1 icarus_stage0_producers.decon2droiEW.wcls_main.params.signal_output_form: "dense" +icarus_stage0_producers.decon2droiEW.wcls_main.logsinks: ["stdout"] +icarus_stage0_producers.decon2droiEW.wcls_main.loglevels: ["debug", "pgraph:info"] +icarus_stage0_producers.decon2droiEW.wcls_main.plugins: [@sequence::icarus_stage0_producers.decon2droiEW.wcls_main.plugins, "WireCellHio", "WireCellPytorch"] +icarus_stage0_producers.decon2droiEW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced.jsonnet"] icarus_stage0_producers.decon2droiWE.wcls_main.inputers: ["wclsRawFrameSource:rfsrc2"] icarus_stage0_producers.decon2droiWE.wcls_main.outputers: ["wclsFrameSaver:spsaver2"] icarus_stage0_producers.decon2droiWE.wcls_main.params.raw_input_label: "daqTPCROI:PHYSCRATEDATATPCWE" icarus_stage0_producers.decon2droiWE.wcls_main.params.tpc_volume_label: 2 icarus_stage0_producers.decon2droiWE.wcls_main.params.signal_output_form: "dense" +icarus_stage0_producers.decon2droiWE.wcls_main.logsinks: ["stdout"] +icarus_stage0_producers.decon2droiWE.wcls_main.loglevels: ["debug", "pgraph:info"] +icarus_stage0_producers.decon2droiWE.wcls_main.plugins: [@sequence::icarus_stage0_producers.decon2droiWE.wcls_main.plugins, "WireCellHio", "WireCellPytorch"] +icarus_stage0_producers.decon2droiWE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced.jsonnet"] icarus_stage0_producers.decon2droiWW.wcls_main.inputers: ["wclsRawFrameSource:rfsrc3"] icarus_stage0_producers.decon2droiWW.wcls_main.outputers: ["wclsFrameSaver:spsaver3"] icarus_stage0_producers.decon2droiWW.wcls_main.params.raw_input_label: "daqTPCROI:PHYSCRATEDATATPCWW" icarus_stage0_producers.decon2droiWW.wcls_main.params.tpc_volume_label: 3 icarus_stage0_producers.decon2droiWW.wcls_main.params.signal_output_form: "dense" +icarus_stage0_producers.decon2droiWW.wcls_main.logsinks: ["stdout"] +icarus_stage0_producers.decon2droiWW.wcls_main.loglevels: ["debug", "pgraph:info"] +icarus_stage0_producers.decon2droiWW.wcls_main.plugins: [@sequence::icarus_stage0_producers.decon2droiWW.wcls_main.plugins, "WireCellHio", "WireCellPytorch"] +icarus_stage0_producers.decon2droiWW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced.jsonnet"] + +# Signal processing values +icarus_stage0_producers.decon2droiEE.wcls_main.structs.gain0: 17.05212 +icarus_stage0_producers.decon2droiEW.wcls_main.structs.gain0: 17.05212 +icarus_stage0_producers.decon2droiWE.wcls_main.structs.gain0: 17.05212 +icarus_stage0_producers.decon2droiWW.wcls_main.structs.gain0: 17.05212 + +icarus_stage0_producers.decon2droiEE.wcls_main.structs.gain1: 12.1420344 +icarus_stage0_producers.decon2droiEW.wcls_main.structs.gain1: 12.1420344 +icarus_stage0_producers.decon2droiWE.wcls_main.structs.gain1: 12.1420344 +icarus_stage0_producers.decon2droiWW.wcls_main.structs.gain1: 12.1420344 + +icarus_stage0_producers.decon2droiEE.wcls_main.structs.gain2: 13.0261362 +icarus_stage0_producers.decon2droiEW.wcls_main.structs.gain2: 13.0261362 +icarus_stage0_producers.decon2droiWE.wcls_main.structs.gain2: 13.0261362 +icarus_stage0_producers.decon2droiWW.wcls_main.structs.gain2: 13.0261362 + + +icarus_stage0_producers.decon2droiEE.wcls_main.structs.shaping0: 1.3 +icarus_stage0_producers.decon2droiEW.wcls_main.structs.shaping0: 1.3 +icarus_stage0_producers.decon2droiWE.wcls_main.structs.shaping0: 1.3 +icarus_stage0_producers.decon2droiWW.wcls_main.structs.shaping0: 1.3 + +icarus_stage0_producers.decon2droiEE.wcls_main.structs.shaping1: 1.45 +icarus_stage0_producers.decon2droiEW.wcls_main.structs.shaping1: 1.45 +icarus_stage0_producers.decon2droiWE.wcls_main.structs.shaping1: 1.45 +icarus_stage0_producers.decon2droiWW.wcls_main.structs.shaping1: 1.45 + +icarus_stage0_producers.decon2droiEE.wcls_main.structs.shaping2: 1.3 +icarus_stage0_producers.decon2droiEW.wcls_main.structs.shaping2: 1.3 +icarus_stage0_producers.decon2droiWE.wcls_main.structs.shaping2: 1.3 +icarus_stage0_producers.decon2droiWW.wcls_main.structs.shaping2: 1.3 + +# Defines final gain, in ADC/e- +icarus_stage0_producers.decon2droiEE.wcls_main.structs.gain_ADC_per_e: @local::gain_ADC_per_e +icarus_stage0_producers.decon2droiEW.wcls_main.structs.gain_ADC_per_e: @local::gain_ADC_per_e +icarus_stage0_producers.decon2droiWE.wcls_main.structs.gain_ADC_per_e: @local::gain_ADC_per_e +icarus_stage0_producers.decon2droiWW.wcls_main.structs.gain_ADC_per_e: @local::gain_ADC_per_e + +# Filter parameters +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Gaus_wide_sigma: @local::Gaus_wide_sigma +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Wiener_tight_U_sigma: @local::Wiener_tight_U_sigma +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Wiener_tight_U_power: @local::Wiener_tight_U_power +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Wiener_tight_V_sigma: @local::Wiener_tight_V_sigma +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Wiener_tight_V_power: @local::Wiener_tight_V_power +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Wiener_tight_W_sigma: @local::Wiener_tight_W_sigma +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Wiener_tight_W_power: @local::Wiener_tight_W_power +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Wire_ind_sigma: @local::Wire_ind_sigma +icarus_stage0_producers.decon2droiEE.wcls_main.structs.Wire_col_sigma: @local::Wire_col_sigma + +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Gaus_wide_sigma: @local::Gaus_wide_sigma +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Wiener_tight_U_sigma: @local::Wiener_tight_U_sigma +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Wiener_tight_U_power: @local::Wiener_tight_U_power +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Wiener_tight_V_sigma: @local::Wiener_tight_V_sigma +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Wiener_tight_V_power: @local::Wiener_tight_V_power +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Wiener_tight_W_sigma: @local::Wiener_tight_W_sigma +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Wiener_tight_W_power: @local::Wiener_tight_W_power +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Wire_ind_sigma: @local::Wire_ind_sigma +icarus_stage0_producers.decon2droiEW.wcls_main.structs.Wire_col_sigma: @local::Wire_col_sigma + +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Gaus_wide_sigma: @local::Gaus_wide_sigma +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Wiener_tight_U_sigma: @local::Wiener_tight_U_sigma +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Wiener_tight_U_power: @local::Wiener_tight_U_power +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Wiener_tight_V_sigma: @local::Wiener_tight_V_sigma +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Wiener_tight_V_power: @local::Wiener_tight_V_power +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Wiener_tight_W_sigma: @local::Wiener_tight_W_sigma +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Wiener_tight_W_power: @local::Wiener_tight_W_power +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Wire_ind_sigma: @local::Wire_ind_sigma +icarus_stage0_producers.decon2droiWE.wcls_main.structs.Wire_col_sigma: @local::Wire_col_sigma + +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Gaus_wide_sigma: @local::Gaus_wide_sigma +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Wiener_tight_U_sigma: @local::Wiener_tight_U_sigma +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Wiener_tight_U_power: @local::Wiener_tight_U_power +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Wiener_tight_V_sigma: @local::Wiener_tight_V_sigma +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Wiener_tight_V_power: @local::Wiener_tight_V_power +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Wiener_tight_W_sigma: @local::Wiener_tight_W_sigma +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Wiener_tight_W_power: @local::Wiener_tight_W_power +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Wire_ind_sigma: @local::Wire_ind_sigma +icarus_stage0_producers.decon2droiWW.wcls_main.structs.Wire_col_sigma: @local::Wire_col_sigma ### Set up to output ROIs from full waveforms icarus_stage0_producers.roifinder1d.WireModuleLabelVec: ["decon1droi:PHYSCRATEDATATPCWW","decon1droi:PHYSCRATEDATATPCWE","decon1droi:PHYSCRATEDATATPCEW","decon1droi:PHYSCRATEDATATPCEE"] @@ -393,10 +501,31 @@ icarus_stage0_producers.gaushit1dTPCWE.CalDataModuleLabel: icarus_stage0_producers.gaushit1dTPCEW.CalDataModuleLabel: "roifinder1d:PHYSCRATEDATATPCEW" icarus_stage0_producers.gaushit1dTPCEE.CalDataModuleLabel: "roifinder1d:PHYSCRATEDATATPCEE" -icarus_stage0_producers.gaushit2dTPCWW.CalDataModuleLabel: "roifinder2d:PHYSCRATEDATATPCWW" -icarus_stage0_producers.gaushit2dTPCWE.CalDataModuleLabel: "roifinder2d:PHYSCRATEDATATPCWE" -icarus_stage0_producers.gaushit2dTPCEW.CalDataModuleLabel: "roifinder2d:PHYSCRATEDATATPCEW" -icarus_stage0_producers.gaushit2dTPCEE.CalDataModuleLabel: "roifinder2d:PHYSCRATEDATATPCEE" +icarus_stage0_producers.gaushit2dTPCWW.CalDataModuleLabel: "decon2droiWW:gauss" +icarus_stage0_producers.gaushit2dTPCWE.CalDataModuleLabel: "decon2droiWE:gauss" +icarus_stage0_producers.gaushit2dTPCEW.CalDataModuleLabel: "decon2droiEW:gauss" +icarus_stage0_producers.gaushit2dTPCEE.CalDataModuleLabel: "decon2droiEE:gauss" + +# Lower thresholds for tighter filter width +icarus_stage0_producers.gaushit2dTPCEE.HitFinderToolVec.CandidateHitsPlane0.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCEE.HitFinderToolVec.CandidateHitsPlane1.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCEE.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCEE.HitFilterAlg.MinPulseHeight: [3., 3., 3.] + +icarus_stage0_producers.gaushit2dTPCEW.HitFinderToolVec.CandidateHitsPlane0.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCEW.HitFinderToolVec.CandidateHitsPlane1.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCEW.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCEW.HitFilterAlg.MinPulseHeight: [3., 3., 3.] + +icarus_stage0_producers.gaushit2dTPCWE.HitFinderToolVec.CandidateHitsPlane0.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCWE.HitFinderToolVec.CandidateHitsPlane1.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCWE.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCWE.HitFilterAlg.MinPulseHeight: [3., 3., 3.] + +icarus_stage0_producers.gaushit2dTPCWW.HitFinderToolVec.CandidateHitsPlane0.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCWW.HitFinderToolVec.CandidateHitsPlane1.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCWW.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 5. +icarus_stage0_producers.gaushit2dTPCWW.HitFilterAlg.MinPulseHeight: [3., 3., 3.] icarus_stage0_producers.gausshitTPCWW.CalDataModuleLabel: "roifinder:PHYSCRATEDATATPCWW" icarus_stage0_producers.gausshitTPCWE.CalDataModuleLabel: "roifinder:PHYSCRATEDATATPCWE" diff --git a/fcl/reco/Stage0/Run2/stage0E_run2_wc_icarus.fcl b/fcl/reco/Stage0/Run2/stage0E_run2_wc_icarus.fcl new file mode 100644 index 000000000..a448be2ca --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0E_run2_wc_icarus.fcl @@ -0,0 +1,10 @@ +#include "stage0_run2_wc_icarus.fcl" + +physics.path: [ + filterdataintegrity, + @sequence::icarus_stage0_PMT, + @sequence::icarus_stage0_CRT, + daqTPCROI, + @sequence::icarus_stage0_multiTPC_2d_TPC_E, + @sequence::icarus_stage0_EastHits2d_TPC +] diff --git a/fcl/reco/Stage0/Run2/stage0E_run2_wc_icarus_mc.fcl b/fcl/reco/Stage0/Run2/stage0E_run2_wc_icarus_mc.fcl new file mode 100644 index 000000000..7c4d7afe4 --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0E_run2_wc_icarus_mc.fcl @@ -0,0 +1,10 @@ +#include "stage0_run2_wc_icarus_mc.fcl" + +# Run the PMT+CRT stage0, initial TPC and East 2D deconv +physics.path: [ @sequence::icarus_stage0_mc_PMT, + MCDecodeTPCROI, + @sequence::icarus_stage0_2d_multiTPC_E, + @sequence::icarus_stage0_mc_crt + ] + +process_name: MCstage0E diff --git a/fcl/reco/Stage0/Run2/stage0E_run2_wcdnn_icarus.fcl b/fcl/reco/Stage0/Run2/stage0E_run2_wcdnn_icarus.fcl new file mode 100644 index 000000000..4dee98c0f --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0E_run2_wcdnn_icarus.fcl @@ -0,0 +1,14 @@ +#include "stage0E_run2_wc_icarus.fcl" +# +# Turn on PCA noise suppression +# physics.producers.MCDecodeTPCROI.DecoderTool.DenoiserType: "pca" + +physics.producers.decon2droiEE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiEW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] + +physics.producers.gaushit2dTPCEE.CalDataModuleLabel: "decon2droiEE:dnnsp" +physics.producers.gaushit2dTPCEW.CalDataModuleLabel: "decon2droiEW:dnnsp" +physics.producers.gaushit2dTPCWE.CalDataModuleLabel: "decon2droiWE:dnnsp" +physics.producers.gaushit2dTPCWW.CalDataModuleLabel: "decon2droiWW:dnnsp" diff --git a/fcl/reco/Stage0/Run2/stage0E_run2_wcdnn_icarus_mc.fcl b/fcl/reco/Stage0/Run2/stage0E_run2_wcdnn_icarus_mc.fcl new file mode 100644 index 000000000..496e730dd --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0E_run2_wcdnn_icarus_mc.fcl @@ -0,0 +1,14 @@ +#include "stage0E_run2_wc_icarus_mc.fcl" +# +# Turn on PCA noise suppression +# physics.producers.MCDecodeTPCROI.DecoderTool.DenoiserType: "pca" + +physics.producers.decon2droiEE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiEW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] + +physics.producers.gaushit2dTPCEE.CalDataModuleLabel: "decon2droiEE:dnnsp" +physics.producers.gaushit2dTPCEW.CalDataModuleLabel: "decon2droiEW:dnnsp" +physics.producers.gaushit2dTPCWE.CalDataModuleLabel: "decon2droiWE:dnnsp" +physics.producers.gaushit2dTPCWW.CalDataModuleLabel: "decon2droiWW:dnnsp" diff --git a/fcl/reco/Stage0/Run2/stage0W_run2_wc_icarus.fcl b/fcl/reco/Stage0/Run2/stage0W_run2_wc_icarus.fcl new file mode 100644 index 000000000..967076819 --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0W_run2_wc_icarus.fcl @@ -0,0 +1,8 @@ +#include "stage0_run2_wc_icarus.fcl" + +physics.path: [ + @sequence::icarus_stage0_multiTPC_2d_TPC_W, + @sequence::icarus_stage0_WestHits2d_TPC +] + +process_name: stage0W diff --git a/fcl/reco/Stage0/Run2/stage0W_run2_wc_icarus_mc.fcl b/fcl/reco/Stage0/Run2/stage0W_run2_wc_icarus_mc.fcl new file mode 100644 index 000000000..5da6bfcc3 --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0W_run2_wc_icarus_mc.fcl @@ -0,0 +1,6 @@ +#include "stage0_run2_wc_icarus_mc.fcl" + +# Run the West 2D deconv and finish TPC stage0 +physics.path: [@sequence::icarus_stage0_2d_multiTPC_W] + +process_name: MCstage0W diff --git a/fcl/reco/Stage0/Run2/stage0W_run2_wcdnn_icarus.fcl b/fcl/reco/Stage0/Run2/stage0W_run2_wcdnn_icarus.fcl new file mode 100644 index 000000000..e1c5866a7 --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0W_run2_wcdnn_icarus.fcl @@ -0,0 +1,14 @@ +#include "stage0W_run2_wc_icarus.fcl" + +# Turn on PCA noise suppression +# physics.producers.MCDecodeTPCROI.DecoderTool.DenoiserType: "pca" + +physics.producers.decon2droiEE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiEW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] + +physics.producers.gaushit2dTPCEE.CalDataModuleLabel: "decon2droiEE:dnnsp" +physics.producers.gaushit2dTPCEW.CalDataModuleLabel: "decon2droiEW:dnnsp" +physics.producers.gaushit2dTPCWE.CalDataModuleLabel: "decon2droiWE:dnnsp" +physics.producers.gaushit2dTPCWW.CalDataModuleLabel: "decon2droiWW:dnnsp" diff --git a/fcl/reco/Stage0/Run2/stage0W_run2_wcdnn_icarus_mc.fcl b/fcl/reco/Stage0/Run2/stage0W_run2_wcdnn_icarus_mc.fcl new file mode 100644 index 000000000..b1661935b --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0W_run2_wcdnn_icarus_mc.fcl @@ -0,0 +1,14 @@ +#include "stage0W_run2_wc_icarus_mc.fcl" + +# Turn on PCA noise suppression +# physics.producers.MCDecodeTPCROI.DecoderTool.DenoiserType: "pca" + +physics.producers.decon2droiEE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiEW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] + +physics.producers.gaushit2dTPCEE.CalDataModuleLabel: "decon2droiEE:dnnsp" +physics.producers.gaushit2dTPCEW.CalDataModuleLabel: "decon2droiEW:dnnsp" +physics.producers.gaushit2dTPCWE.CalDataModuleLabel: "decon2droiWE:dnnsp" +physics.producers.gaushit2dTPCWW.CalDataModuleLabel: "decon2droiWW:dnnsp" diff --git a/fcl/reco/Stage0/Run2/stage0_run2_wc_icarus.fcl b/fcl/reco/Stage0/Run2/stage0_run2_wc_icarus.fcl index 50fc99937..e9bc6bc56 100644 --- a/fcl/reco/Stage0/Run2/stage0_run2_wc_icarus.fcl +++ b/fcl/reco/Stage0/Run2/stage0_run2_wc_icarus.fcl @@ -32,39 +32,6 @@ outputs.rootOutput.outputCommands: ["keep *_*_*_*", "keep *_daq_ICARUSTriggerV*_*" ] -# Add below as per Gray Putnam (should we move to jsonnet files?) -physics.producers.decon2droiEE.wcls_main.structs.gain0: 17.05212 -physics.producers.decon2droiEW.wcls_main.structs.gain0: 17.05212 -physics.producers.decon2droiWE.wcls_main.structs.gain0: 17.05212 -physics.producers.decon2droiWW.wcls_main.structs.gain0: 17.05212 - - -physics.producers.decon2droiEE.wcls_main.structs.gain1: 12.1420344 -physics.producers.decon2droiEW.wcls_main.structs.gain1: 12.1420344 -physics.producers.decon2droiWE.wcls_main.structs.gain1: 12.1420344 -physics.producers.decon2droiWW.wcls_main.structs.gain1: 12.1420344 - -physics.producers.decon2droiEE.wcls_main.structs.gain2: 13.0261362 -physics.producers.decon2droiEW.wcls_main.structs.gain2: 13.0261362 -physics.producers.decon2droiWE.wcls_main.structs.gain2: 13.0261362 -physics.producers.decon2droiWW.wcls_main.structs.gain2: 13.0261362 - - -physics.producers.decon2droiEE.wcls_main.structs.shaping0: 1.3 -physics.producers.decon2droiEW.wcls_main.structs.shaping0: 1.3 -physics.producers.decon2droiWE.wcls_main.structs.shaping0: 1.3 -physics.producers.decon2droiWW.wcls_main.structs.shaping0: 1.3 - -physics.producers.decon2droiEE.wcls_main.structs.shaping1: 1.45 -physics.producers.decon2droiEW.wcls_main.structs.shaping1: 1.45 -physics.producers.decon2droiWE.wcls_main.structs.shaping1: 1.45 -physics.producers.decon2droiWW.wcls_main.structs.shaping1: 1.45 - -physics.producers.decon2droiEE.wcls_main.structs.shaping2: 1.3 -physics.producers.decon2droiEW.wcls_main.structs.shaping2: 1.3 -physics.producers.decon2droiWE.wcls_main.structs.shaping2: 1.3 -physics.producers.decon2droiWW.wcls_main.structs.shaping2: 1.3 - services.message.destinations : { STDCOUT: diff --git a/fcl/reco/Stage0/Run2/stage0_run2_wc_icarus_mc.fcl b/fcl/reco/Stage0/Run2/stage0_run2_wc_icarus_mc.fcl index 71f48966e..26fef30e6 100644 --- a/fcl/reco/Stage0/Run2/stage0_run2_wc_icarus_mc.fcl +++ b/fcl/reco/Stage0/Run2/stage0_run2_wc_icarus_mc.fcl @@ -40,45 +40,10 @@ physics.producers.ophit.InputModule: "opdaq" physics.producers.MCDecodeTPCROI.FragmentsLabelVec: ["daq:TPCWW","daq:TPCWE","daq:TPCEW","daq:TPCEE"] physics.producers.MCDecodeTPCROI.OutInstanceLabelVec: ["PHYSCRATEDATATPCWW", "PHYSCRATEDATATPCWE", "PHYSCRATEDATATPCEW", "PHYSCRATEDATATPCEE"] -physics.producers.decon1droi.RawDigitLabelVec: ["MCDecodeTPCROI:PHYSCRATEDATATPCWW","MCDecodeTPCROI:PHYSCRATEDATATPCWE","MCDecodeTPCROI:PHYSCRATEDATATPCEW","MCDecodeTPCROI:PHYSCRATEDATATPCEE"] - physics.producers.decon2droiEE.wcls_main.params.raw_input_label: "MCDecodeTPCROI:PHYSCRATEDATATPCEE" physics.producers.decon2droiEW.wcls_main.params.raw_input_label: "MCDecodeTPCROI:PHYSCRATEDATATPCEW" physics.producers.decon2droiWE.wcls_main.params.raw_input_label: "MCDecodeTPCROI:PHYSCRATEDATATPCWE" physics.producers.decon2droiWW.wcls_main.params.raw_input_label: "MCDecodeTPCROI:PHYSCRATEDATATPCWW" -# As per Gray Putname... -physics.producers.decon2droiEE.wcls_main.structs.gain0: 17.05212 -physics.producers.decon2droiEW.wcls_main.structs.gain0: 17.05212 -physics.producers.decon2droiWE.wcls_main.structs.gain0: 17.05212 -physics.producers.decon2droiWW.wcls_main.structs.gain0: 17.05212 - - -physics.producers.decon2droiEE.wcls_main.structs.gain1: 12.1420344 -physics.producers.decon2droiEW.wcls_main.structs.gain1: 12.1420344 -physics.producers.decon2droiWE.wcls_main.structs.gain1: 12.1420344 -physics.producers.decon2droiWW.wcls_main.structs.gain1: 12.1420344 - -physics.producers.decon2droiEE.wcls_main.structs.gain2: 13.0261362 -physics.producers.decon2droiEW.wcls_main.structs.gain2: 13.0261362 -physics.producers.decon2droiWE.wcls_main.structs.gain2: 13.0261362 -physics.producers.decon2droiWW.wcls_main.structs.gain2: 13.0261362 - - -physics.producers.decon2droiEE.wcls_main.structs.shaping0: 1.3 -physics.producers.decon2droiEW.wcls_main.structs.shaping0: 1.3 -physics.producers.decon2droiWE.wcls_main.structs.shaping0: 1.3 -physics.producers.decon2droiWW.wcls_main.structs.shaping0: 1.3 - -physics.producers.decon2droiEE.wcls_main.structs.shaping1: 1.45 -physics.producers.decon2droiEW.wcls_main.structs.shaping1: 1.45 -physics.producers.decon2droiWE.wcls_main.structs.shaping1: 1.45 -physics.producers.decon2droiWW.wcls_main.structs.shaping1: 1.45 - -physics.producers.decon2droiEE.wcls_main.structs.shaping2: 1.3 -physics.producers.decon2droiEW.wcls_main.structs.shaping2: 1.3 -physics.producers.decon2droiWE.wcls_main.structs.shaping2: 1.3 -physics.producers.decon2droiWW.wcls_main.structs.shaping2: 1.3 - # restore legacy G4 labels physics.producers.mcophit.SimPhotonsProducer: "largeant" diff --git a/fcl/reco/Stage0/Run2/stage0_run2_wcdnn_icarus.fcl b/fcl/reco/Stage0/Run2/stage0_run2_wcdnn_icarus.fcl new file mode 100644 index 000000000..70b694516 --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0_run2_wcdnn_icarus.fcl @@ -0,0 +1,14 @@ +#include "stage0_run2_wc_icarus.fcl" + +# Turn on PCA noise suppression +# physics.producers.MCDecodeTPCROI.DecoderTool.DenoiserType: "pca" + +physics.producers.decon2droiEE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiEW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] + +physics.producers.gaushit2dTPCEE.CalDataModuleLabel: "decon2droiEE:dnnsp" +physics.producers.gaushit2dTPCEW.CalDataModuleLabel: "decon2droiEW:dnnsp" +physics.producers.gaushit2dTPCWE.CalDataModuleLabel: "decon2droiWE:dnnsp" +physics.producers.gaushit2dTPCWW.CalDataModuleLabel: "decon2droiWW:dnnsp" diff --git a/fcl/reco/Stage0/Run2/stage0_run2_wcdnn_icarus_mc.fcl b/fcl/reco/Stage0/Run2/stage0_run2_wcdnn_icarus_mc.fcl new file mode 100644 index 000000000..cdfb878da --- /dev/null +++ b/fcl/reco/Stage0/Run2/stage0_run2_wcdnn_icarus_mc.fcl @@ -0,0 +1,40 @@ +#include "stage0_run2_wc_icarus_mc.fcl" + +# Turn on PCA noise suppression +# physics.producers.MCDecodeTPCROI.DecoderTool.DenoiserType: "pca" + +physics.producers.decon2droiEE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiEW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWE.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] +physics.producers.decon2droiWW.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet"] + +physics.producers.gaushit2dTPCEE.CalDataModuleLabel: "decon2droiEE:dnnsp" +physics.producers.gaushit2dTPCEW.CalDataModuleLabel: "decon2droiEW:dnnsp" +physics.producers.gaushit2dTPCWE.CalDataModuleLabel: "decon2droiWE:dnnsp" +physics.producers.gaushit2dTPCWW.CalDataModuleLabel: "decon2droiWW:dnnsp" + +physics.producers.gaushit2dTPCEE.HitFinderToolVec.CandidateHitsPlane0.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCEE.HitFinderToolVec.CandidateHitsPlane1.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCEE.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCEE.HitFilterAlg.MinPulseHeight: [2., 2., 2.] + +physics.producers.gaushit2dTPCEW.HitFinderToolVec.CandidateHitsPlane0.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCEW.HitFinderToolVec.CandidateHitsPlane1.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCEW.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCEW.HitFilterAlg.MinPulseHeight: [2., 2., 2.] + +physics.producers.gaushit2dTPCWE.HitFinderToolVec.CandidateHitsPlane0.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCWE.HitFinderToolVec.CandidateHitsPlane1.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCWE.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCWE.HitFilterAlg.MinPulseHeight: [2., 2., 2.] + +physics.producers.gaushit2dTPCWW.HitFinderToolVec.CandidateHitsPlane0.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCWW.HitFinderToolVec.CandidateHitsPlane1.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCWW.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 2.5 +physics.producers.gaushit2dTPCWW.HitFilterAlg.MinPulseHeight: [2., 2., 2.] + +physics.analyzers.spanaE.HitProducers: ["gaushit2dTPCEE", "gaushit2dTPCEW"] +physics.analyzers.spanaE.WireProducers: ["decon2droiEE:dnnsp", "decon2droiEW:dnnsp"] + +physics.analyzers.spanaW.HitProducers: ["gaushit2dTPCWE", "gaushit2dTPCWW"] +physics.analyzers.spanaW.WireProducers: ["decon2droiWE:dnnsp", "decon2droiWW:dnnsp"] diff --git a/icaruscode/TPC/CMakeLists.txt b/icaruscode/TPC/CMakeLists.txt index 20d12e1bc..3bca0be5a 100644 --- a/icaruscode/TPC/CMakeLists.txt +++ b/icaruscode/TPC/CMakeLists.txt @@ -8,3 +8,4 @@ add_subdirectory(SignalProcessing) add_subdirectory(Simulation) add_subdirectory(Tracking) add_subdirectory(ICARUSWireCell) +add_subdirectory(SPAna) diff --git a/icaruscode/TPC/ICARUSWireCell/detsimmodules_wirecell_ICARUS.fcl b/icaruscode/TPC/ICARUSWireCell/detsimmodules_wirecell_ICARUS.fcl index fb6d9d486..8debc90bc 100644 --- a/icaruscode/TPC/ICARUSWireCell/detsimmodules_wirecell_ICARUS.fcl +++ b/icaruscode/TPC/ICARUSWireCell/detsimmodules_wirecell_ICARUS.fcl @@ -12,9 +12,10 @@ icarus_simwire_wirecell: { tool_type: WCLS apps: ["Pgrapher"] - // logsinks: ["stdout"] + logsinks: ["stdout"] // loglevels: ["magnify:debug"] - plugins: ["WireCellPgraph", "WireCellGen","WireCellSio","WireCellRoot","WireCellLarsoft"] + loglevels: ["debug", "pgraph:info"] + plugins: ["WireCellPgraph", "WireCellGen","WireCellSio","WireCellRoot","WireCellLarsoft", "WireCellHio"] // needs to be found via your WIRECELL_PATH configs: ["pgrapher/experiment/icarus/wcls-multitpc-sim-drift-simchannel-refactored.jsonnet"] // Contract note: these exact "type:name" must be used to identify @@ -75,6 +76,13 @@ icarus_simwire_wirecell: } } +# SimChannel only (SConly) CONFIG +icarus_simwire_wirecell_SConly: @local::icarus_simwire_wirecell +icarus_simwire_wirecell_SConly.wcls_main.configs: ["pgrapher/experiment/icarus/wcls-multitpc-sim-drift-simchannel-true.jsonnet"] +icarus_simwire_wirecell_SConly.wcls_main.outputers:[ + "wclsDepoFluxWriter:postdrift" +] + # TUNED FR CONFIG icarus_simwire_wirecell_fitSR: @local::icarus_simwire_wirecell # Add in the ER tail diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/chndb-base.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/chndb-base.jsonnet index 02ba66873..48f34af7a 100644 --- a/icaruscode/TPC/ICARUSWireCell/icarus/chndb-base.jsonnet +++ b/icaruscode/TPC/ICARUSWireCell/icarus/chndb-base.jsonnet @@ -28,7 +28,7 @@ function(params, anode, field, n, rms_cuts=[]) // Externally determined "bad" channels. - bad: [], + bad: [167, 1795, 1796, 1824, 1826, 1828, 1830, 1860, 1888, 1889, 1890, 1891, 1922, 2013, 2806, 3725, 3743, 4262, 5572, 8523, 8524, 8525, 8526, 8527, 8945, 9180, 9387, 9391, 9664, 10138, 12051, 14028, 14403, 15307, 15737, 16017, 16018, 16019, 16020, 16021, 16022, 16023, 16024, 16025, 16026, 16027, 16028, 16029, 16030, 16031, 16937, 17201, 19580, 23225, 23457, 23984, 23985, 25662, 27219, 27231, 27335, 27967, 28640, 28804, 28805, 28806, 28807, 28808, 28809, 28810, 30496, 30497, 30498, 30499, 30500, 30501, 30502, 30503, 30504, 30505, 30506, 30507, 30508, 30509, 30510, 30511, 30512, 30513, 30514, 30515, 30516, 30517, 30518, 30519, 30520, 30521, 30522, 30523, 30524, 30525, 30526, 30527, 30705, 31316, 31317, 31726, 34220, 38294, 41501, 41975, 42851, 43231, 43518, 44051, 44410, 49216, 49792, 49793, 49794, 49795, 49796, 49797, 49798, 49799, 49800, 49801, 49802, 49803, 49804, 49805, 49806, 49807, 49808, 49809, 49810, 49811, 49812, 49813, 49814, 49815, 49816, 49817, 49818, 49819, 49820, 49821, 49822, 49823, 50080, 50081, 50082, 50083, 50084, 50085, 50086, 50087, 50088, 50089, 50090, 50091, 50092, 50093, 50094, 50095, 50096, 50097, 50098, 50099, 50100, 50101, 50102, 50103, 50104, 50105, 50106, 50107, 50108, 50109, 50110, 50111, 50673, 51900], // Overide defaults for specific channels. If an info is // mentioned for a particular channel in multiple objects in this @@ -41,7 +41,7 @@ function(params, anode, field, n, rms_cuts=[]) // repeat values found here in subsequent entries unless you // wish to change them. { - channels: util.anode_channels(n), + channels: util.anode_channels_twofaced(n), nominal_baseline: 2048.0, // adc count gain_correction: 1.0, // unitless response_offset: 0.0, // ticks? diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/dnnroi.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/dnnroi.jsonnet new file mode 100644 index 000000000..173e8ecc9 --- /dev/null +++ b/icaruscode/TPC/ICARUSWireCell/icarus/dnnroi.jsonnet @@ -0,0 +1,111 @@ +// This produces a function to configure DNN-ROI for one APA given +// anode and torch service (ts) objects. +// +// The prefix is prepended to all internal node names if uniqueness +// beyond anode ID is needed. The output_scale allows for an ad-hoc +// scaling of dnnroi output. The U and W planes will go through +// dnnroi while hte W plane will be shunted. What comes out will be a +// unified frame with frame tag "dnnspN" where "N" is the anode ID. + +local wc = import "wirecell.jsonnet"; +local pg = import "pgraph.jsonnet"; + +function (anode, ts_u, ts_v, prefix="dnnroi", output_scale=1, nchunk_u=2, nchunk_v=4) + local apaid = anode.data.ident; + local prename = prefix + std.toString(apaid); + local intags = ['loose_lf%d'%apaid, 'mp2_roi%d'%apaid, + 'mp3_roi%d'%apaid]; + + local dnnroi_u = pg.pnode({ + type: "DNNROIFinding", + name: prename+"u", + data: { + anode: wc.tn(anode), + plane: 0, + nticks: 4096, + intags: intags, + decon_charge_tag: "decon%d" %apaid, + outtag: "dnnsp%du"%apaid, + output_scale: output_scale, + input_scale: 0.00025, + forward: wc.tn(ts_u), + tick_per_slice: 4, + nchunks: nchunk_u + } + }, nin=1, nout=1, uses=[ts_u, anode]); + local dnnroi_v = pg.pnode({ + type: "DNNROIFinding", + name: prename+"v", + data: { + anode: wc.tn(anode), + plane: 1, + nticks: 4096, + intags: intags, + decon_charge_tag: "decon%d" %apaid, + outtag: "dnnsp%dv"%apaid, + output_scale: output_scale, + input_scale: 0.00025, + forward: wc.tn(ts_v), + tick_per_slice: 4, + nchunks: nchunk_v + } + }, nin=1, nout=1, uses=[ts_v, anode]); + local dnnroi_w = pg.pnode({ + type: "PlaneSelector", + name: prename+"w", + data: { + anode: wc.tn(anode), + plane: 2, + tags: ["gauss%d"%apaid], + tag_rules: [{ + frame: {".*":"DNNROIFinding"}, + trace: {["gauss%d"%apaid]:"dnnsp%dw"%apaid}, + }], + } + }, nin=1, nout=1, uses=[anode]); + + local dnnpipes = [dnnroi_u, dnnroi_v, dnnroi_w]; + local dnnfanout = pg.pnode({ + type: "FrameFanout", + name: prename, + data: { + multiplicity: 3 + } + }, nin=1, nout=3); + + local dnnfanin = pg.pnode({ + type: "FrameFanin", + name: prename, + data: { + multiplicity: 3, + tag_rules: [{ + frame: {".*": "dnnsp%d%s" % [apaid,plane]}, + trace: {".*": "dnnsp%d%s" % [apaid,plane]}, + } for plane in ["u", "v", "w"]] + }, + }, nin=3, nout=1); + + local retagger = pg.pnode({ + type: "Retagger", + name: 'dnnroi%d' % apaid, + data: { + // Note: retagger keeps tag_rules an array to be like frame fanin/fanout. + tag_rules: [{ + // Retagger also handles "frame" and "trace" like fanin/fanout + // merge separately all traces like gaussN to gauss. + frame: { + ".*": "dnnsp%d" % apaid + }, + merge: { + ["dnnsp%d." % apaid]: "dnnsp%d" % apaid + }, + }], + }, + }, nin=1, nout=1); + + pg.intern(innodes=[dnnfanout], + outnodes=[retagger], + centernodes=dnnpipes+[dnnfanin], + edges=[pg.edge(dnnfanout, dnnpipes[ind], ind, 0) for ind in [0,1,2]] + + [pg.edge(dnnpipes[ind], dnnfanin, 0, ind) for ind in [0,1,2]] + + [pg.edge(dnnfanin, retagger, 0, 0)]) diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/funcs.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/funcs.jsonnet index 5f2f0be20..d0e6fee8e 100644 --- a/icaruscode/TPC/ICARUSWireCell/icarus/funcs.jsonnet +++ b/icaruscode/TPC/ICARUSWireCell/icarus/funcs.jsonnet @@ -25,6 +25,12 @@ local g = import 'pgraph.jsonnet'; anode_channels(n):: std.flattenArrays([std.range(startch[n][w], startch[n][w]+wireplanes[w]) for w in std.range(0,2)]), // anode_channels(n):: std.range(1056 * (n % 2) + 13312 * (n - n % 2) / 2, 1056 * (n % 2 + 1) - 1 + 13312 * (n - n % 2) / 2) + std.range(1056 * 2 + 13312 * (n - n % 2) / 2, 13312 - 1 + 13312 * (n - n % 2) / 2), + // Channels on each anode, in the "two-faced" configuration + anode_channels_twofaced(n):: std.range(startch[2*n][0], startch[2*n][0]+wireplanes[0]) + + std.range(startch[2*n+1][0], startch[2*n+1][0]+wireplanes[0]) + + std.range(startch[2*n][1], startch[2*n][1]+wireplanes[1]) + + std.range(startch[2*n][2], startch[2*n][2]+wireplanes[2]), + // Return the number of split (1 or 2) for an anode anode_split(ident):: (ident%100 - ident%10)/10, diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/nf.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/nf.jsonnet index cc8b758ce..8cecf585a 100644 --- a/icaruscode/TPC/ICARUSWireCell/icarus/nf.jsonnet +++ b/icaruscode/TPC/ICARUSWireCell/icarus/nf.jsonnet @@ -3,34 +3,8 @@ local g = import 'pgraph.jsonnet'; local wc = import 'wirecell.jsonnet'; -local default_dft = { type: 'FftwDFT' }; - -function(params, anode, chndbobj, tools, name='', dft=default_dft) +function(anode, chndbobj, tools, name='') { - - local single = { - type: 'icarusOneChannelNoise', - name: name, - uses: [dft, chndbobj, anode, tools.rc_resp], - data: { - noisedb: wc.tn(chndbobj), - anode: wc.tn(anode), - rcresp: wc.tn(tools.rc_resp), - dft: wc.tn(dft), - }, - }, - local grouped = { - type: 'mbCoherentNoiseSub', - name: name, - uses: [dft, chndbobj, anode], - data: { - noisedb: wc.tn(chndbobj), - anode: wc.tn(anode), - dft: wc.tn(dft), - rms_threshold: 0.0, - }, - }, - local obnf = g.pnode({ type: 'OmnibusNoiseFilter', name: name, @@ -43,18 +17,16 @@ function(params, anode, chndbobj, tools, name='', dft=default_dft) // only when the channelmask is merged to `bad` // maskmap: {sticky: "bad", ledge: "bad", noisy: "bad"}, channel_filters: [ - wc.tn(single), ], grouped_filters: [ - // wc.tn(grouped), ], channel_status_filters: [ ], noisedb: wc.tn(chndbobj), - intraces: 'orig%d' % anode.data.ident, // frame tag get all traces - outtraces: 'raw%d' % anode.data.ident, + intraces: 'orig%d' % anode.data.ident, + outtraces: 'orig%d' % anode.data.ident, }, - }, uses=[chndbobj, anode, tools.rc_resp, single, grouped], nin=1, nout=1), + }, uses=[chndbobj, anode], nin=1, nout=1), pipe: g.pipeline([obnf], name=name), diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/params_twofaced.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/params_twofaced.jsonnet new file mode 100644 index 000000000..e31a81f37 --- /dev/null +++ b/icaruscode/TPC/ICARUSWireCell/icarus/params_twofaced.jsonnet @@ -0,0 +1,142 @@ +// This file specifies the paramter configuration for the ICARUS detector. It +// inherit from the base params.jsonnet and override the relevant parameters + +local wc = import "wirecell.jsonnet"; +local base = import "pgrapher/common/params.jsonnet"; + +base { + det : { + + // define the 4 APAs. This must use the coordinate system + // defined by the wire geometry file. + // + // The "faces" is consumed by, at least, the Drifter and + // AnodePlane. The "wires" number is used to set + // AnodePlane.ident used to lookup the anode in WireSchema. + // It corresponds to the anode number. + // + // Also see: + // lar -c dump_icarus_geometry.fcl + // wirecell-util wire-volumes icarus-wires-dualanode.json.bz2 + // to help with defining these parameters. + + + // The "response" plane is where the field response functions + // start. Garfield calcualtions start somewhere relative to + // something, here's where that is made concrete. This MUST + // match what field response functions also used. + response_plane: 10*wc.cm, // relative to collection wires + + + // Each wire gets identified with the id number in the geometry + // file. The horizontal induction is split into two, both + // sides are enumerated by "s", while "n" counts the physical anodes. + // NOTE: the actual physical volumes in ICARUS are only 4 + + local xanode = [-359.33*wc.cm, -61.1*wc.cm, 61.1*wc.cm, 359.33*wc.cm], + local offset_response = [if a%2==0 then +10*wc.cm else -10*wc.cm for a in std.range(0,3)], + local xresponse = [xanode[a] + offset_response[a] for a in std.range(0,3)], + local xcathode = [-210.29*wc.cm, -210.29*wc.cm, 210.29*wc.cm, 210.29*wc.cm], + volumes : [ + { + local anode = a, // physical anode number + + wires: (anode), + name: "anode%d"%(anode), + + faces: [ + { + anode: xanode[a], + response: xresponse[a], + cathode: xcathode[a], + }, + { + anode: xanode[a], + response: xresponse[a], + cathode: xcathode[a], + } + ], + } for a in std.range(0,3) + ], + + // This describes some rough, overall bounding box. It's not + // directly needed but can be useful on the Jsonnet side, for + // example when defining some simple kinematics. It is + // represented by a ray going from extreme corners of a + // rectangular solid. Again "wirecell-util wires-info" helps + // to choose something. //FIXME -- ARE CORRECT? + bounds : { + tail: wc.point(-3.65, -1.7, -9.1, wc.m), + head: wc.point(+3.65, +1.4, +8.8, wc.m), + } + }, + + daq: super.daq { + tick: 0.4*wc.us, // 2.5 MHz + nticks: 4096, + }, + + adc: super.adc { + // fix baseline at 2048 (induction), 400 (collection) + baselines: [1650.0*wc.millivolt, 1650.0*wc.millivolt, 322.3*wc.millivolt], + + // From ICARUS paper: https://iopscience.iop.org/article/10.1088/1748-0221/13/12/P12007/pdf + //check (values taken from the FE calibration shown in pg. 7 of the paper) + fullscale: [0.8*wc.millivolt, 3.3*wc.volt], + }, + + elec: [super.elec { + type: "WarmElecResponse", + // Old values: + // ICARUS nominal: 17.8075*wc.mV/wc.fC // 0.027 fC/(ADC*us) + // Match data ADC values (docdb 25161): 14.9654*wc.mV/wc.fC, // 0.0321 fC/(ADC*us) + gain: 14.9654*wc.mV/wc.fC, // 0.0321 fC/(ADC*us) + shaping: 1.3*wc.us, + postgain: 1.0, + start: 0, + }, for _ in [0, 1, 2]], + + + sim: super.sim { + + // For running in LArSoft, the simulation must be in fixed time mode. + fixed: true, + + // The "absolute" time (ie, in G4 time) that the lower edge of + // of final readout tick #0 should correspond to. This is a + // "fixed" notion. + local tick0_time = -340*wc.us, // TriggerOffsetTPC from detectorclocks_icarus.fcl + + // Open the ductor's gate a bit early. + local response_time_offset = $.det.response_plane / $.lar.drift_speed, + local response_nticks = wc.roundToInt(response_time_offset / $.daq.tick), + + ductor : { + nticks: $.daq.nticks + response_nticks, + readout_time: self.nticks * $.daq.tick, + start_time: tick0_time - response_time_offset, + }, + + // To counter the enlarged duration of the ductor, a Reframer + // chops off the little early, extra time. Note, tags depend on how + reframer: { + tbin: response_nticks, + nticks: $.daq.nticks, + } + + }, + + files: { + wires: "icarus_wire_twofaced.json.bz2", + + fields: ["icarus_fnal_fit_ks_P0nom.json.bz2"], + + // noise: ["icarus_noise_model_int_TPCEE.json.bz2","icarus_noise_model_int_TPCEW.json.bz2","icarus_noise_model_int_TPCWE.json.bz2","icarus_noise_model_int_TPCWW.json.bz2"], + // coherent_noise: ["icarus_noise_model_coh_TPCEE.json.bz2","icarus_noise_model_coh_TPCEW.json.bz2","icarus_noise_model_coh_TPCWE.json.bz2","icarus_noise_model_coh_TPCWW.json.bz2"], + wiregroups: "icarus_group_to_channel_map.json.bz2", + noisegroups: ["icarus_noise_model_int_by_board_TPCEE.json.bz2","icarus_noise_model_int_by_board_TPCEW.json.bz2","icarus_noise_model_int_by_board_TPCWE.json.bz2","icarus_noise_model_int_by_board_TPCWW.json.bz2","icarus_noise_model_coh_by_board_TPCEE.json.bz2","icarus_noise_model_coh_by_board_TPCEW.json.bz2","icarus_noise_model_coh_by_board_TPCWE.json.bz2","icarus_noise_model_coh_by_board_TPCWW.json.bz2"], + chresp: null, + }, + +} + diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/sp-filters.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/sp-filters.jsonnet index 695b2ec4e..607773a5d 100644 --- a/icaruscode/TPC/ICARUSWireCell/icarus/sp-filters.jsonnet +++ b/icaruscode/TPC/ICARUSWireCell/icarus/sp-filters.jsonnet @@ -43,19 +43,22 @@ local wf(name, data={}) = { lf('ROI_loose_lf', { tau: 0.0025 * wc.megahertz }), // 0.0025 hf('Gaus_tight'), - hf('Gaus_wide', { sigma: 0.12 * wc.megahertz }), + hf('Gaus_wide', { + sigma: std.extVar("Gaus_wide_sigma")*wc.megahertz //0.12 * wc.megahertz + }), hf('Wiener_tight_U', { - sigma: 0.148788 * wc.megahertz, - power: 3.76194, + sigma: std.extVar('Wiener_tight_U_sigma')*wc.megahertz, // 0.148788 * wc.megahertz, + power: std.extVar('Wiener_tight_U_power'), //3.76194, }), hf("Wiener_tight_V", { - sigma: 0.1596568 * wc.megahertz, - power: 4.36125 }), + sigma: std.extVar('Wiener_tight_V_sigma')*wc.megahertz, // 0.1596568 * wc.megahertz, + power: std.extVar('Wiener_tight_V_power'), // 4.36125 + }), hf('Wiener_tight_W', { - sigma: 0.13623 * wc.megahertz, - power: 3.35324, + sigma: std.extVar('Wiener_tight_W_sigma')*wc.megahertz, // 0.13623 * wc.megahertz, + power: std.extVar('Wiener_tight_W_power'), // 3.35324, }), hf('Wiener_wide_U', { @@ -71,6 +74,10 @@ local wf(name, data={}) = { power: 4.37928, }), - wf('Wire_ind', { sigma: 1.0 / wc.sqrtpi * 0.75 }), - wf('Wire_col', { sigma: 1.0 / wc.sqrtpi * 3.0 }), + wf('Wire_ind', { + sigma: 1.0 / wc.sqrtpi * std.extVar("Wire_ind_sigma") // 0.75 + }), + wf('Wire_col', { + sigma: 1.0 / wc.sqrtpi * std.extVar("Wire_col_sigma") // 3.0 + }), ] diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/sp.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/sp.jsonnet index ace6dd8be..10650e6bd 100644 --- a/icaruscode/TPC/ICARUSWireCell/icarus/sp.jsonnet +++ b/icaruscode/TPC/ICARUSWireCell/icarus/sp.jsonnet @@ -30,7 +30,7 @@ function(params, tools, override = {}) { elecresponse : wc.tn(tools.elec_resp[2]), postgain: 1, // default 1.2 ADC_mV: 4096 / (3300.0 * wc.mV), // default 4096/2000 - troi_col_th_factor: 5.0, // default 5 + troi_col_th_factor: 3.0, // default 5 troi_ind_th_factor: 3.0, // default 3 lroi_rebin: 6, // default 6 lroi_th_factor: 3.5, // default 3.5 diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/splat.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/splat.jsonnet new file mode 100644 index 000000000..14e05416b --- /dev/null +++ b/icaruscode/TPC/ICARUSWireCell/icarus/splat.jsonnet @@ -0,0 +1,62 @@ +local wc = import "wirecell.jsonnet"; +local g = import 'pgraph.jsonnet'; + +// The approximated sim+sigproc +function(params, tools, anode, name=null) { + local sufix = name, + local sp = g.pnode({ + type: 'DepoFluxSplat', + name: sufix, + data: { + anode: wc.tn(anode), + field_response: wc.tn(tools.field), // for speed and origin + sparse: true, + tick: params.daq.tick, + window_start: -340 * wc.us, // TriggerOffsetTPC from detectorclocks_icarus.fcl + window_duration: params.daq.readout_time, + reference_time: -1500 * wc.us - self.window_start, // G4RefTime from detectorclocks_icarus.fcl less window start as per Brett Viren + // Run wirecell-gen morse-* to find these numbers that match the extra + // spread the sigproc induces. + "smear_long": [ + 6.6, + 6.6, + 6.6, + ], + "smear_tran": [ + 1.4, + 1.4, + 0.4, + ] + }, + }, nin=1, nout=1, uses=[anode, tools.field]), + local hio = g.pnode({ + type: 'HDF5FrameTap', + name: 'hio_true%s' % name, + data: { + anode: wc.tn(anode), + trace_tags: ['deposplat%s' % name], + filename: "wc-true-%s.h5" % name, + chunk: [0, 0], // ncol, nrow + gzip: 2, + high_throughput: true, + }, + }, nin=1, nout=1), + local rt = g.pnode({ + type: 'Retagger', + name: sufix, + data: { + // Note: retagger keeps tag_rules an array to be like frame fanin/fanout. + tag_rules: [{ + // Retagger also handles "frame" and "trace" like fanin/fanout + // merge separately all traces like gaussN to gauss. + frame: { + ".*": "deposplat%s" % sufix + }, + merge: { + ".*": "deposplat%s" % sufix + }, + }], + }, + }, nin=1, nout=1), + ret: g.pipeline([sp, rt, hio], name=sp.name), +}.ret diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet new file mode 100644 index 000000000..923cdeed4 --- /dev/null +++ b/icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig-twofaced-dnnroi.jsonnet @@ -0,0 +1,390 @@ +// This is a main entry point to configure a WC/LS job that applies +// noise filtering and signal processing to existing RawDigits. The +// FHiCL is expected to provide the following parameters as attributes +// in the "params" structure. +// +// epoch: the hardware noise fix expoch: "before", "after", "dynamic" or "perfect" +// reality: whether we are running on "data" or "sim"ulation. +// raw_input_label: the art::Event inputTag for the input RawDigit +// +// see the .fcl of the same name for an example Version 2 +// +// Manual testing, eg: +// +// jsonnet -V reality=data -V epoch=dynamic -V raw_input_label=daq \\ +// -V signal_output_form=sparse \\ +// -J cfg cfg/pgrapher/experiment/uboone/wcls-nf-sp.jsonnet +// +// jsonnet -V reality=sim -V epoch=perfect -V raw_input_label=daq \\ +// -V signal_output_form=sparse \\ +// -J cfg cfg/pgrapher/experiment/uboone/wcls-nf-sp.jsonnet + + +local epoch = std.extVar('epoch'); // eg "dynamic", "after", "before", "perfect" +local reality = std.extVar('reality'); +local sigoutform = std.extVar('signal_output_form'); // eg "sparse" or "dense" + +local wc = import 'wirecell.jsonnet'; +local g = import 'pgraph.jsonnet'; + +local raw_input_label = std.extVar('raw_input_label'); // eg "daq" +local volume_label = std.extVar('tpc_volume_label'); // eg "",0,1,2,3 +local volume = if volume_label == '' then -1 else std.parseInt(volume_label); + +// local data_params = import 'params.jsonnet'; +// local simu_params = import 'simparams.jsonnet'; +// local params_init = if reality == 'data' then data_params else simu_params; +local params_twofaced = import 'pgrapher/experiment/icarus/params_twofaced.jsonnet'; + +# Load the sim-params, overwrite the volume config for the two-faced version +local base_params = import 'pgrapher/experiment/icarus/simparams.jsonnet'; +local base = base_params + params_twofaced; + +// load the electronics response parameters +local er_params = [ + { + gain: std.extVar('gain0')*wc.mV/wc.fC, + shaping: std.extVar('shaping0')*wc.us, + }, + + { + gain: std.extVar('gain1')*wc.mV/wc.fC, + shaping: std.extVar('shaping1')*wc.us, + }, + + { + gain: std.extVar('gain2')*wc.mV/wc.fC, + shaping: std.extVar('shaping2')*wc.us, + }, +]; + + +local params = base { + files: super.files { + fields: [ std.extVar('files_fields'), ], + chresp: null, + }, + + rc_resp: if std.extVar('file_rcresp') != "" then + { + // "icarus_fnal_rc_tail.json" + filename: std.extVar('file_rcresp'), + postgain: 1.0, + start: 0.0, + tick: 0.4*wc.us, + nticks: params.daq.nticks,// 4255, + type: "JsonElecResponse", + } + else super.rc_resp, + + elec: std.mapWithIndex(function (n, eparam) + super.elec[n] + { + gain: eparam.gain, + shaping: eparam.shaping, + }, er_params), + +}; + +// local tools_maker = import 'pgrapher/common/tools.jsonnet'; +local tools_maker = import 'pgrapher/experiment/icarus/icarus_tools.jsonnet'; +local tools = tools_maker(params); + + +local wcls_maker = import 'pgrapher/ui/wcls/nodes.jsonnet'; +local wcls = wcls_maker(params, tools); + +local sp_maker = import 'pgrapher/experiment/icarus/sp.jsonnet'; + +//local chndbm = chndb_maker(params, tools); +//local chndb = if epoch == "dynamic" then chndbm.wcls_multi(name="") else chndbm.wct(epoch); + + +// Collect the WC/LS input converters for use below. Make sure the +// "name" argument matches what is used in the FHiCL that loads this +// file. In particular if there is no ":" in the inputer then name +// must be the emtpy string. +local wcls_input = { + adc_digits: g.pnode({ + type: 'wclsRawFrameSource', + name: 'rfsrc%d' %volume, // to use multiple wirecell instances in a fhicl job + data: { + art_tag: raw_input_label, + frame_tags: ['orig'], // this is a WCT designator + tick: params.daq.tick, + // nticks: params.daq.nticks, + }, + }, nin=0, nout=1), + +}; + +// Collect all the wc/ls output converters for use below. Note the +// "name" MUST match what is used in theh "outputers" parameter in the +// FHiCL that loads this file. + +local this_anode = tools.anodes[volume]; + +local wcls_output = { + // The noise filtered "ADC" values. These are truncated for + // art::Event but left as floats for the WCT SP. Note, the tag + // "raw" is somewhat historical as the output is not equivalent to + // "raw data". + nf_digits: g.pnode({ + type: 'wclsFrameSaver', + name: 'nfsaver', + data: { + // anode: wc.tn(tools.anode), + anode: wc.tn(this_anode), + digitize: true, // true means save as RawDigit, else recob::Wire + frame_tags: ['raw'], + // nticks: params.daq.nticks, + chanmaskmaps: ['bad'], + plane_map: { + "1": 3, // front induction: WireCell::kULayer -> geo::kH (1 -> 3) + "2": 1, // middle induction: WireCell:kVLayer -> geo::kV (2 -> 1) + "4": 0, // collection: WireCell::kWLayer -> geo::kU (4 -> 0) + }, + }, + }, nin=1, nout=1, uses=[this_anode]), + + + // The output of signal processing. Note, there are two signal + // sets each created with its own filter. The "gauss" one is best + // for charge reconstruction, the "wiener" is best for S/N + // separation. Both are used in downstream WC code. + sp_signals: g.pnode({ + type: 'wclsFrameSaver', + name: 'spsaver%d' %volume, // to use multiple wirecell instances in a fhicl job + data: { + // anode: wc.tn(tools.anode), + anode: wc.tn(this_anode), + digitize: false, // true means save as RawDigit, else recob::Wire + // frame_tags: ['gauss', 'wiener', 'looseLf','shrinkROI','extendROI'], + // frame_scale: [0.1, 0.1, 0.1], + // frame_tags: ['gauss','wiener','looseLf','shrinkROI','extendROI','mp3ROI','mp2ROI', 'cleanupROI'], + // frame_scale: [0.009,0.009,0.009,0.009,0.009,0.009,0.009,0.009], + + frame_tags: ['dnnsp'], + frame_scale: [std.extVar('gain_ADC_per_e')], + + // nticks: params.daq.nticks, + chanmaskmaps: ['bad'], + nticks: -1, + plane_map: { + "1": 3, // front induction: WireCell::kULayer -> geo::kH (1 -> 3) + "2": 1, // middle induction: WireCell:kVLayer -> geo::kV (2 -> 1) + "4": 0, // collection: WireCell::kWLayer -> geo::kU (4 -> 0) + }, + }, + }, nin=1, nout=1, uses=[this_anode]), + + h5io: g.pnode({ + type: 'HDF5FrameTap', + name: 'hio_sp%d' % volume, + data: { + anode: wc.tn(this_anode), + trace_tags: ['gauss' + , 'wiener' + , 'tightLf' + , 'looseLf' + , 'decon' + , 'cleanupROI' + , 'breakROI1' + , 'breakROI2' + , 'shrinkROI' + , 'extendROI' + , 'mp3ROI' + , 'mp2ROI' + , 'dnnsp' + ], + filename: "wc-sp-%d.h5" % volume, + chunk: [0, 0], // ncol, nrow + gzip: 0, + high_throughput: true, + }, + }, nin=1, nout=1, uses=[this_anode]), + +}; + +// local perfect = import 'chndb-perfect.jsonnet'; +local base = import 'pgrapher/experiment/icarus/chndb-base.jsonnet'; +local chndb = [{ + type: 'OmniChannelNoiseDB', + name: 'ocndbperfect%d' % n, + // data: perfect(params, tools.anodes[n], tools.field, n), + data: base(params, tools.anodes[n], tools.field, n){dft:wc.tn(tools.dft)}, + uses: [tools.anodes[n], tools.field, tools.dft], // pnode extension +} for n in std.range(0, std.length(tools.anodes) - 1)]; + +local nf_maker = import 'pgrapher/experiment/icarus/nf.jsonnet'; +local nf_pipes = [nf_maker(tools.anodes[n], chndb[n], tools, name='nf%d' % n) for n in std.range(0, std.length(tools.anodes) - 1)]; + +local sp_override = { // assume all tages sets in base sp.jsonnet + sparse: sigoutform == 'sparse', + use_roi_refinement: true, + use_roi_debug_mode: true, + save_negtive_charge: true, + wiener_tag: "", + // gauss_tag: "", + tight_lf_tag: "", + // loose_lf_tag: "", + break_roi_loop1_tag: "", + break_roi_loop2_tag: "", + shrink_roi_tag: "", + extend_roi_tag: "", + // decon_charge_tag: "", + cleanup_roi_tag: "", + // mp2_roi_tag: "", + // mp3_roi_tag: "", + use_multi_plane_protection: true, + mp_tick_resolution: 8, + process_planes: [0, 1, 2], + isWrapped: true, + nwires_separate_planes: [ + [1056, 1056], [5600], [5600] + ], + troi_col_th_factor: std.parseJson(std.extVar('col_threshold_factor'))*1.0, // multiply by 1 to make into float + troi_ind_th_factor: std.parseJson(std.extVar('ind_threshold_factor'))*1.0 +}; +local sp = sp_maker(params, tools, sp_override); +local sp_pipes = [sp.make_sigproc(a) for a in tools.anodes]; + +local util = import 'pgrapher/experiment/icarus/funcs.jsonnet'; +local chsel_pipes = [ + g.pnode({ + type: 'ChannelSelector', + name: 'chsel%d' % n, + data: { + channels: util.anode_channels_twofaced(n), + //tags: ['orig%d' % n], // traces tag + }, + }, nin=1, nout=1) + for n in std.range(0, std.length(tools.anodes) - 1) +]; + +local magoutput = 'icarus-data-check.root'; +local magnify = import 'pgrapher/experiment/icarus/magnify-sinks.jsonnet'; +local magnifyio = magnify(tools, magoutput); + +local dnnroi = import 'pgrapher/experiment/icarus/dnnroi.jsonnet'; +local ts_u = { + type: "TorchService", + name: "dnnroi_u", + data: { + model: "NNsv5/plane0.ts", + device: "cpu", + concurrency: 1, + }, +}; + +local ts_v = { + type: "TorchService", + name: "dnnroi_v", + data: { + model: "NNsv5/plane1.ts", + device: "cpu", + concurrency: 1, + }, +}; + +local nfsp_pipes = [ + g.pipeline([ + chsel_pipes[n], + // magnifyio.orig_pipe[n], + + nf_pipes[n], + // magnifyio.raw_pipe[n], + sp_pipes[n], + dnnroi(tools.anodes[n], ts_u, ts_v), + // magnifyio.decon_pipe[n], + // magnifyio.threshold_pipe[n], + // magnifyio.debug_pipe[n], // use_roi_debug_mode: true in sp.jsonnet + ], + 'nfsp_pipe_%d' % n) + for n in [volume] +]; + +local fanout_tag_rules = [ + { + frame: { + '.*': 'orig%d' % tools.anodes[n].data.ident, + }, + trace: { + // fake doing Nmult SP pipelines + //orig: ['wiener', 'gauss'], + //'.*': 'orig', + }, + } + for n in [volume] + ]; + +local fanin_tag_rules = [ + { + frame: { + //['number%d' % n]: ['output%d' % n, 'output'], + '.*': 'framefanin', + }, + trace: { + ['extend_roi%d'%ind]:'extend_roi%d'%ind, + ['shrink_roi%d'%ind]:'shrink_roi%d'%ind, + // ['break_roi_2nd%d'%ind]:'break_roi_2nd%d'%ind, + // ['break_roi_1st%d'%ind]:'break_roi_1st%d'%ind, + ['cleanup_roi%d'%ind]:'cleanup_roi%d'%ind, + ['mp2_roi%d'%ind]:'mp2_roi%d'%ind, + ['mp3_roi%d'%ind]:'mp3_roi%d'%ind, + ['gauss%d'%ind]:'gauss%d'%ind, + 'dnnsp\\d':'dnnsp%d'%ind, + ['wiener%d'%ind]:'wiener%d'%ind, + // ['threshold%d'%ind]:'threshold%d'%ind, + // ['tight_lf%d'%ind]:'tight_lf%d'%ind, + ['loose_lf%d'%ind]:'loose_lf%d'%ind, + // ['decon%d'%ind]:'decon%d'%ind, + }, + + } + for ind in [this_anode.data.ident] + ]; +local fanpipe = util.fanpipe('FrameFanout', nfsp_pipes, 'FrameFanin', 'nfsp%d' % volume, [], fanout_tag_rules, fanin_tag_rules); + +local retagger = g.pnode({ + type: 'Retagger', + data: { + // Note: retagger keeps tag_rules an array to be like frame fanin/fanout. + tag_rules: [{ + // Retagger also handles "frame" and "trace" like fanin/fanout + // merge separately all traces like gaussXYZ to gauss. + frame: { + '.*': 'retagger', + }, + merge: { + 'dnnsp\\d': 'dnnsp', + 'gauss\\d': 'gauss', + 'wiener\\d': 'wiener', + // 'tight_lf\\d': 'tightLf', + 'loose_lf\\d': 'looseLf', + // 'decon\\d': 'decon', + 'cleanup_roi\\d': 'cleanupROI', + // 'break_roi_1st\\d': 'breakROI1', + // 'break_roi_2nd\\d': 'breakROI2', + 'shrink_roi\\d': 'shrinkROI', + 'extend_roi\\d': 'extendROI', + 'mp3_roi\\d': 'mp3ROI', + 'mp2_roi\\d': 'mp2ROI', + }, + }], + }, +}, nin=1, nout=1); + +local sink = g.pnode({ type: 'DumpFrames' }, nin=1, nout=0); + +// local graph = g.pipeline([wcls_input.adc_digits, chsel_pipes[volume], sp_pipes[volume], dnnroi(this_anode, ts_u, ts_v, output_scale=1), retagger, wcls_output.h5io, wcls_output.sp_signals, sink]); +local graph = g.pipeline([wcls_input.adc_digits, fanpipe, retagger, wcls_output.sp_signals, sink]); + +local app = { + type: 'Pgrapher', + data: { + edges: g.edges(graph), + }, +}; + +// Finally, the configuration sequence +g.uses(graph) + [app] diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig-twofaced.jsonnet similarity index 69% rename from icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig.jsonnet rename to icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig-twofaced.jsonnet index 5355d1779..a804569da 100644 --- a/icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig.jsonnet +++ b/icaruscode/TPC/ICARUSWireCell/icarus/wcls-decode-to-sig-twofaced.jsonnet @@ -34,7 +34,11 @@ local volume = if volume_label == '' then -1 else std.parseInt(volume_label); // local data_params = import 'params.jsonnet'; // local simu_params = import 'simparams.jsonnet'; // local params_init = if reality == 'data' then data_params else simu_params; -local base = import 'pgrapher/experiment/icarus/simparams.jsonnet'; +local params_twofaced = import 'pgrapher/experiment/icarus/params_twofaced.jsonnet'; + +# Load the sim-params, overwrite the volume config for the two-faced version +local base_params = import 'pgrapher/experiment/icarus/simparams.jsonnet'; +local base = base_params + params_twofaced; // load the electronics response parameters local er_params = [ @@ -116,14 +120,9 @@ local wcls_input = { // Collect all the wc/ls output converters for use below. Note the // "name" MUST match what is used in theh "outputers" parameter in the // FHiCL that loads this file. -local mega_anode = { - type: 'MegaAnodePlane', - name: 'meganodes', - data: { - anodes_tn: if volume != -1 then [wc.tn(a) for a in tools.anodes[2*volume:2*(volume+1)]] // single volume - else [wc.tn(anode) for anode in tools.anodes], // all volumes - }, -}; + +local this_anode = tools.anodes[volume]; + local wcls_output = { // The noise filtered "ADC" values. These are truncated for // art::Event but left as floats for the WCT SP. Note, the tag @@ -134,13 +133,18 @@ local wcls_output = { name: 'nfsaver', data: { // anode: wc.tn(tools.anode), - anode: wc.tn(mega_anode), + anode: wc.tn(this_anode), digitize: true, // true means save as RawDigit, else recob::Wire frame_tags: ['raw'], // nticks: params.daq.nticks, chanmaskmaps: ['bad'], + plane_map: { + "1": 3, // front induction: WireCell::kULayer -> geo::kH (1 -> 3) + "2": 1, // middle induction: WireCell:kVLayer -> geo::kV (2 -> 1) + "4": 0, // collection: WireCell::kWLayer -> geo::kU (4 -> 0) + }, }, - }, nin=1, nout=1, uses=[mega_anode]), + }, nin=1, nout=1, uses=[this_anode]), // The output of signal processing. Note, there are two signal @@ -152,17 +156,52 @@ local wcls_output = { name: 'spsaver%d' %volume, // to use multiple wirecell instances in a fhicl job data: { // anode: wc.tn(tools.anode), - anode: wc.tn(mega_anode), + anode: wc.tn(this_anode), digitize: false, // true means save as RawDigit, else recob::Wire - // frame_tags: ['gauss', 'wiener', 'looseLf'], + // frame_tags: ['gauss', 'wiener', 'looseLf','shrinkROI','extendROI'], // frame_scale: [0.1, 0.1, 0.1], - frame_tags: ['decon','looseLf'], - frame_scale: [0.009,0.009], + // frame_tags: ['gauss','wiener','looseLf','shrinkROI','extendROI','mp3ROI','mp2ROI', 'cleanupROI'], + // frame_scale: [0.009,0.009,0.009,0.009,0.009,0.009,0.009,0.009], + + frame_tags: ['gauss'], + frame_scale: [std.extVar('gain_ADC_per_e')], + // nticks: params.daq.nticks, - chanmaskmaps: [], + chanmaskmaps: ['bad'], nticks: -1, + plane_map: { + "1": 3, // front induction: WireCell::kULayer -> geo::kH (1 -> 3) + "2": 1, // middle induction: WireCell:kVLayer -> geo::kV (2 -> 1) + "4": 0, // collection: WireCell::kWLayer -> geo::kU (4 -> 0) + }, }, - }, nin=1, nout=1, uses=[mega_anode]), + }, nin=1, nout=1, uses=[this_anode]), + + h5io: g.pnode({ + type: 'HDF5FrameTap', + name: 'hio_sp%d' % volume, + data: { + anode: wc.tn(this_anode), + trace_tags: ['gauss' + , 'wiener' + , 'tightLf' + , 'looseLf' + , 'decon' + , 'cleanupROI' + , 'breakROI1' + , 'breakROI2' + , 'shrinkROI' + , 'extendROI' + , 'mp3ROI' + , 'mp2ROI' + ], + filename: "wc-sp-%d.h5" % volume, + chunk: [0, 0], // ncol, nrow + gzip: 0, + high_throughput: true, + }, + }, nin=1, nout=1, uses=[this_anode]), + }; // local perfect = import 'chndb-perfect.jsonnet'; @@ -176,24 +215,34 @@ local chndb = [{ } for n in std.range(0, std.length(tools.anodes) - 1)]; local nf_maker = import 'pgrapher/experiment/icarus/nf.jsonnet'; -local nf_pipes = [nf_maker(params, tools.anodes[n], chndb[n], tools, name='nf%d' % n) for n in std.range(0, std.length(tools.anodes) - 1)]; +local nf_pipes = [nf_maker(tools.anodes[n], chndb[n], tools, name='nf%d' % n) for n in std.range(0, std.length(tools.anodes) - 1)]; local sp_override = { // assume all tages sets in base sp.jsonnet sparse: sigoutform == 'sparse', + use_roi_refinement: true, + use_roi_debug_mode: false, + save_negtive_charge: true, // wiener_tag: "", // gauss_tag: "", - use_roi_refinement: false, - use_roi_debug_mode: true, tight_lf_tag: "", - // loose_lf_tag: "", - cleanup_roi_tag: "", + loose_lf_tag: "", break_roi_loop1_tag: "", break_roi_loop2_tag: "", shrink_roi_tag: "", extend_roi_tag: "", - // m_decon_charge_tag: "", + decon_charge_tag: "", + cleanup_roi_tag: "", + mp2_roi_tag: "", + mp3_roi_tag: "", use_multi_plane_protection: false, - mp_tick_resolution: 10, + mp_tick_resolution: 8, + process_planes: [0, 1, 2], + isWrapped: true, + nwires_separate_planes: [ + [1056, 1056], [5600], [5600] + ], + troi_col_th_factor: std.parseJson(std.extVar('col_threshold_factor'))*1.0, // multiply by 1 to make into float + troi_ind_th_factor: std.parseJson(std.extVar('ind_threshold_factor'))*1.0 }; local sp = sp_maker(params, tools, sp_override); local sp_pipes = [sp.make_sigproc(a) for a in tools.anodes]; @@ -204,7 +253,7 @@ local chsel_pipes = [ type: 'ChannelSelector', name: 'chsel%d' % n, data: { - channels: util.anode_channels(n), + channels: util.anode_channels_twofaced(n), //tags: ['orig%d' % n], // traces tag }, }, nin=1, nout=1) @@ -220,7 +269,7 @@ local nfsp_pipes = [ chsel_pipes[n], // magnifyio.orig_pipe[n], - // nf_pipes[n], + nf_pipes[n], // magnifyio.raw_pipe[n], sp_pipes[n], @@ -254,18 +303,25 @@ local fanin_tag_rules = [ '.*': 'framefanin', }, trace: { + ['extend_roi%d'%ind]:'extend_roi%d'%ind, + ['shrink_roi%d'%ind]:'shrink_roi%d'%ind, + // ['break_roi_2nd%d'%ind]:'break_roi_2nd%d'%ind, + // ['break_roi_1st%d'%ind]:'break_roi_1st%d'%ind, + ['cleanup_roi%d'%ind]:'cleanup_roi%d'%ind, + ['mp2_roi%d'%ind]:'mp2_roi%d'%ind, + ['mp3_roi%d'%ind]:'mp3_roi%d'%ind, ['gauss%d'%ind]:'gauss%d'%ind, ['wiener%d'%ind]:'wiener%d'%ind, - ['threshold%d'%ind]:'threshold%d'%ind, - // ['tight_lf%d'%ind]:'tight_lf%d'%ind, + // ['threshold%d'%ind]:'threshold%d'%ind, + // ['tight_lf%d'%ind]:'tight_lf%d'%ind, ['loose_lf%d'%ind]:'loose_lf%d'%ind, - ['decon%d'%ind]:'decon%d'%ind, + // ['decon%d'%ind]:'decon%d'%ind, }, } for ind in anode_ident ]; -local fanpipe = util.fanpipe('FrameFanout', nfsp_pipes, 'FrameFanin', 'nfsp', [], fanout_tag_rules, fanin_tag_rules); +local fanpipe = util.fanpipe('FrameFanout', nfsp_pipes, 'FrameFanin', 'nfsp%d' % volume, [], fanout_tag_rules, fanin_tag_rules); local retagger = g.pnode({ type: 'Retagger', @@ -278,11 +334,18 @@ local retagger = g.pnode({ '.*': 'retagger', }, merge: { - 'gauss\\d\\d\\d': 'gauss', - 'wiener\\d\\d\\d': 'wiener', - // 'tight_lf\\d\\d\\d': 'tightLf', - 'loose_lf\\d\\d\\d': 'looseLf', - 'decon\\d\\d\\d': 'decon', + 'gauss\\d': 'gauss', + 'wiener\\d': 'wiener', + // 'tight_lf\\d': 'tightLf', + 'loose_lf\\d': 'looseLf', + // 'decon\\d': 'decon', + 'cleanup_roi\\d': 'cleanupROI', + // 'break_roi_1st\\d': 'breakROI1', + // 'break_roi_2nd\\d': 'breakROI2', + 'shrink_roi\\d': 'shrinkROI', + 'extend_roi\\d': 'extendROI', + 'mp3_roi\\d': 'mp3ROI', + 'mp2_roi\\d': 'mp2ROI', }, }], }, @@ -290,7 +353,6 @@ local retagger = g.pnode({ local sink = g.pnode({ type: 'DumpFrames' }, nin=1, nout=0); - local graph = g.pipeline([wcls_input.adc_digits, fanpipe, retagger, wcls_output.sp_signals, sink]); local app = { diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/wcls-multitpc-sim-drift-simchannel-true.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/wcls-multitpc-sim-drift-simchannel-true.jsonnet new file mode 100644 index 000000000..35775075d --- /dev/null +++ b/icaruscode/TPC/ICARUSWireCell/icarus/wcls-multitpc-sim-drift-simchannel-true.jsonnet @@ -0,0 +1,178 @@ +local g = import 'pgraph.jsonnet'; +local f = import 'pgrapher/common/funcs.jsonnet'; +local wc = import 'wirecell.jsonnet'; + +local io = import 'pgrapher/common/fileio.jsonnet'; +local tools_maker = import 'pgrapher/experiment/icarus/icarus_tools.jsonnet'; +local base = import 'pgrapher/experiment/icarus/simparams.jsonnet'; +local splat = import 'pgrapher/experiment/icarus/splat.jsonnet'; + +// load the electronics response parameters +local er_params = [ + { + gain: std.extVar('gain0')*wc.mV/wc.fC, + shaping: std.extVar('shaping0')*wc.us, + }, + + { + gain: std.extVar('gain1')*wc.mV/wc.fC, + shaping: std.extVar('shaping1')*wc.us, + }, + + { + gain: std.extVar('gain2')*wc.mV/wc.fC, + shaping: std.extVar('shaping2')*wc.us, + }, +]; + +local params = base { + lar: super.lar { + // Longitudinal diffusion constant + DL: std.extVar('DL') * wc.cm2 / wc.ns, + // Transverse diffusion constant + DT: std.extVar('DT') * wc.cm2 / wc.ns, + // Electron lifetime + lifetime: std.extVar('lifetime') * wc.us, + // Electron drift speed, assumes a certain applied E-field + // drift_speed: std.extVar('driftSpeed') * wc.mm / wc.us, + }, + files: super.files { + fields: [ std.extVar('files_fields'), ], + }, + + rc_resp: if std.extVar('file_rcresp') != "" then + { + // "icarus_fnal_rc_tail.json" + filename: std.extVar('file_rcresp'), + postgain: 1.0, + start: 0.0, + tick: 0.4*wc.us, + nticks: 4255, + type: "JsonElecResponse", + rc_layers: 1 + } + else super.rc_resp, + + elec: std.mapWithIndex(function (n, eparam) + super.elec[n] + { + gain: eparam.gain, + shaping: eparam.shaping, + }, er_params), +}; + +local tools = tools_maker(params); + +local sim_maker = import 'pgrapher/experiment/icarus/sim.jsonnet'; +local sim = sim_maker(params, tools); + +local nanodes = std.length(tools.anodes); +local anode_iota = std.range(0, nanodes - 1); + +local output = 'wct-sim-ideal-sig.npz'; + + +local wcls_maker = import 'pgrapher/ui/wcls/nodes.jsonnet'; +local wcls = wcls_maker(params, tools); + +//Haiwang DepoSetSource Implementation: +local wcls_input = { + depos: wcls.input.depos(name="", art_tag="IonAndScint"), + deposet: g.pnode({ + type: 'wclsSimDepoSetSource', + name: "electron", + data: { + model: "", + scale: -1, //scale is -1 to correct a sign error in the SimDepoSource converter. + art_tag: "ionization", //name of upstream art producer of depos "label:instance:processName" + assn_art_tag: "", + id_is_track: false, // Use this for "id-is-index" in the output + }, + }, nin=0, nout=1), +}; + +// Collect all the wc/ls output converters for use below. Note the +// "name" MUST match what is used in theh "outputers" parameter in the +// FHiCL that loads this file. +local mega_anode = { + type: 'MegaAnodePlane', + name: 'meganodes', + data: { + anodes_tn: [wc.tn(anode) for anode in tools.anodes], + }, +}; + +// A ``duo'' anode consists of two ``splits'' +local duoanodes = [ +{ + type: 'MegaAnodePlane', + name: 'duoanode%d' %n, + data: { + // anodes_tn: ["AnodePlane:anode110", "AnodePlane:anode120"], + anodes_tn: [wc.tn(a) for a in tools.anodes[2*n:2*(n+1)]], + // anodes_tn: [wc.tn(tools.anodes[2*n]), wc.tn(tools.anodes[2*n+1])], + }, +} +for n in std.range(0,3)]; +local volname = ["EE", "EW", "WE", "WW"]; + +local anode_names = ["EES", "EEN", "EWS", "EWN", "WES", "WEN", "WWS", "WWN"]; + +local drifter = sim.drifter; +local setdrifter = g.pnode({ + type: 'DepoSetDrifter', + data: { + drifter: "Drifter" + } + }, nin=1, nout=1, + uses=[drifter]); + +local wcls_simchannel_sink = + g.pnode({ + type: 'wclsDepoFluxWriter', + name: 'postdrift', + data: { + anodes: [wc.tn(anode) for anode in tools.anodes], + field_response: wc.tn(tools.field), + + // time binning + tick: params.daq.tick, + window_start: -340 * wc.us, // TriggerOffsetTPC from detectorclocks_icarus.fcl + window_duration: params.daq.readout_time, + + nsigma: 3.0, + + reference_time: -1500 * wc.us - self.window_start, // G4RefTime from detectorclocks_icarus.fcl less window start as per Brett Viren + + smear_long: 0.0, + smear_tran: 0.0, + + time_offsets: [std.extVar('time_offset_u') * wc.us, std.extVar('time_offset_v') * wc.us, std.extVar('time_offset_y') * wc.us], + + // input from art::Event + sed_label: 'largeant:TPCActive', + + // output to art::Event + simchan_label: 'simpleSC', + }, + }, nin=1, nout=1, uses=tools.anodes+[tools.field]); + +local util = import 'pgrapher/experiment/icarus/funcs.jsonnet'; + +local deposplats = [splat(params, tools, tools.anodes[n], name="A%d" % n) for n in anode_iota] ; +local makesplat = util.fanpipe("DepoSetFanout", deposplats, "FrameFanin"); + +local sink = sim.frame_sink; + +local graph = g.pipeline([wcls_input.deposet, setdrifter, wcls_simchannel_sink, makesplat, sink]); + +local app = { + type: 'Pgrapher', + data: { + edges: g.edges(graph), + }, +}; + + +// Finally, the configuration sequence which is emitted. + +g.uses(graph) + [app] diff --git a/icaruscode/TPC/ICARUSWireCell/icarus/wcls-multitpc-sim-drift-simchannel.jsonnet b/icaruscode/TPC/ICARUSWireCell/icarus/wcls-multitpc-sim-drift-simchannel.jsonnet index b5516b88f..6a85eb784 100644 --- a/icaruscode/TPC/ICARUSWireCell/icarus/wcls-multitpc-sim-drift-simchannel.jsonnet +++ b/icaruscode/TPC/ICARUSWireCell/icarus/wcls-multitpc-sim-drift-simchannel.jsonnet @@ -150,6 +150,35 @@ local wcls_output = { }, nin=1, nout=1, uses=[duoanodes[n]]) for n in std.range(0,3)], + // Output to H5 files + h5io: [g.pnode({ + type: 'HDF5FrameTap', + name: 'hio_sp%d' % n, + data: { + anode: wc.tn(duoanodes[n]), + trace_tags: ['TPC%s' %volname[n] + , 'loose_lf%s' % volname[n] + , 'tight_lf%s' % volname[n] + , 'cleanup_roi%s' % volname[n] + , 'break_roi_1st%s' % volname[n] + , 'break_roi_2nd%s' % volname[n] + , 'shrink_roi%s' % volname[n] + , 'extend_roi%s' % volname[n] + , 'mp3_roi%s' % volname[n] + , 'mp2_roi%s' % volname[n] + , 'decon_charge%s' % volname[n] + , 'wiener%s' % volname[n] + , 'gauss%s' % volname[n]], + filename: "wc-rd-%s.h5" % volname[n], + chunk: [0, 0], // ncol, nrow + gzip: 2, + high_throughput: true, + }, + }, nin=1, nout=1), + for n in std.range(0, 3)], + + + // The noise filtered "ADC" values. These are truncated for // art::Event but left as floats for the WCT SP. Note, the tag // "raw" is somewhat historical as the output is not equivalent to @@ -196,10 +225,6 @@ local chndb = [{ // local nf_pipes = [nf_maker(params, tools.anodes[n], chndb_pipes[n]) for n in std.range(0, std.length(tools.anodes)-1)]; // local nf_pipes = [nf_maker(params, tools.anodes[n], chndb[n], n, name='nf%d' % n) for n in anode_iota]; -local sp_maker = import 'pgrapher/experiment/icarus/sp.jsonnet'; -local sp = sp_maker(params, tools); -local sp_pipes = [sp.make_sigproc(a) for a in tools.anodes]; - local rng = tools.random; local wcls_simchannel_sink_old = g.pnode({ @@ -329,8 +354,12 @@ local frame_summers = [ }, }, nin=2, nout=1) for n in std.range(0, 3)]; -local actpipes = [g.pipeline([noises[n], coh_noises[n], digitizers[n], /*retaggers[n],*/ wcls_output.sim_digits[n]], name="noise-digitizer%d" %n) for n in std.range(0,3)]; local util = import 'pgrapher/experiment/icarus/funcs.jsonnet'; + +// local actpipes = [g.pipeline([noises[n], coh_noises[n], digitizers[n], /*retaggers[n],*/ wcls_output.sim_digits[n]], name="noise-digitizer%d" %n) for n in std.range(0,3)]; + +local actpipes = [g.pipeline([noises[n], coh_noises[n], digitizers[n], wcls_output.h5io[n], wcls_output.sim_digits[n]], name="noise-digitizer%d" %n) for n in std.range(0,3)]; + local outtags = ['orig%d' % n for n in std.range(0, 3)]; local pipe_reducer = util.fansummer('DepoSetFanout', analog_pipes, frame_summers, actpipes, 'FrameFanin', 'fansummer', outtags); @@ -356,7 +385,7 @@ local sink = sim.frame_sink; // local graph = g.pipeline([wcls_input.depos, drifter, wcls_simchannel_sink.simchannels, bagger, pipe_reducer, retagger, wcls_output.sim_digits, sink]); //local graph = g.pipeline([wcls_input.depos, drifter, wcls_simchannel_sink, bagger, pipe_reducer, sink]); -local graph = g.pipeline([wcls_input.deposet, setdrifter, wcls_simchannel_sink_old, wcls_simchannel_sink, pipe_reducer, sink]); +local graph = g.pipeline([wcls_input.deposet, setdrifter, wcls_simchannel_sink, wcls_simchannel_sink_old, pipe_reducer, sink]); local app = { type: 'Pgrapher', diff --git a/icaruscode/TPC/ICARUSWireCell/wcls-decode-to-sig-base.fcl b/icaruscode/TPC/ICARUSWireCell/wcls-decode-to-sig-base.fcl index 0c8cfaa85..c97b8fa12 100644 --- a/icaruscode/TPC/ICARUSWireCell/wcls-decode-to-sig-base.fcl +++ b/icaruscode/TPC/ICARUSWireCell/wcls-decode-to-sig-base.fcl @@ -54,6 +54,10 @@ standard_wirecell_sigproc: signal_output_form: "sparse" file_rcresp: "icarus_fnal_rc_tail.json" + + // ROI threshold, defined as factor times RMS noise. Per ind/col. + col_threshold_factor: 3.0 + ind_threshold_factor: 3.0 } } } diff --git a/icaruscode/TPC/SPAna/CMakeLists.txt b/icaruscode/TPC/SPAna/CMakeLists.txt new file mode 100644 index 000000000..4bfa8badd --- /dev/null +++ b/icaruscode/TPC/SPAna/CMakeLists.txt @@ -0,0 +1,34 @@ +cet_build_plugin(SPAna art::module + LIBRARIES + art::Framework_Core + art::Framework_Services_Registry + art_root_io::TFileService_service + art::Framework_Principal + art::Persistency_Common + art::Utilities canvas::canvas + cetlib::cetlib cetlib_except::cetlib_except + ROOT::X3d + art_root_io::tfile_support + art_root_io::art_root_io + art_root_io::dict + lardataobj::AnalysisBase + lardataobj::RecoBase + larcorealg::Geometry + larcore::Geometry_Geometry_service + larcorealg::GeoAlgo + sbnobj::Common_Reco + larcorealg::GeoAlgo + lardata::DetectorInfoServices_DetectorPropertiesServiceStandard_service + lardataalg::DetectorInfo + ROOT::Minuit + sbnobj::Common_Calibration_dict + larevt::SpaceCharge + nug4::ParticleNavigation + larsim::MCCheater_BackTrackerService_service + larsim::MCCheater_ParticleInventoryService_service + larsim::Simulation +) + +install_headers() +install_source() +install_fhicl() diff --git a/icaruscode/TPC/SPAna/SPAna_module.cc b/icaruscode/TPC/SPAna/SPAna_module.cc new file mode 100644 index 000000000..4a8330e49 --- /dev/null +++ b/icaruscode/TPC/SPAna/SPAna_module.cc @@ -0,0 +1,682 @@ +//////////////////////////////////////////////////////////////////////// +// Class: SPAna +// Plugin Type: analyzer (Unknown Unknown) +// File: SPAna_module.cc +// +// Generated at Sun Jun 23 05:48:04 2024 by Gray Putnam using cetskelgen +// from version . +//////////////////////////////////////////////////////////////////////// + +#include "art/Framework/Core/EDAnalyzer.h" +#include "art/Framework/Core/ModuleMacros.h" +#include "art/Framework/Principal/Event.h" +#include "art/Framework/Principal/Handle.h" +#include "art/Framework/Principal/Run.h" +#include "art/Framework/Principal/SubRun.h" +#include "canvas/Utilities/InputTag.h" +#include "fhiclcpp/ParameterSet.h" +#include "messagefacility/MessageLogger/MessageLogger.h" + +#include "TTree.h" +#include "art_root_io/TFileService.h" + +#include "canvas/Persistency/Common/FindMany.h" +#include "canvas/Persistency/Common/FindManyP.h" +#include "canvas/Persistency/Common/FindOne.h" +#include "canvas/Persistency/Common/FindOneP.h" +#include "canvas/Persistency/Common/Ptr.h" +#include "canvas/Persistency/Common/PtrVector.h" + +#include "lardataobj/RecoBase/SpacePoint.h" +#include "lardataobj/RecoBase/Wire.h" +#include "lardataobj/RecoBase/Hit.h" +#include "lardataobj/RecoBase/Cluster.h" +#include "lardataobj/AnalysisBase/T0.h" +#include "nusimdata/SimulationBase/MCParticle.h" +#include "nusimdata/SimulationBase/MCTruth.h" + +#include "sbnobj/Common/Calibration/TrackCaloSkimmerObj.h" + +#include "larevt/SpaceCharge/SpaceCharge.h" +#include "larevt/SpaceChargeServices/SpaceChargeService.h" +#include "lardataalg/DetectorInfo/DetectorPropertiesStandard.h" +#include "lardata/DetectorInfoServices/DetectorPropertiesService.h" +#include "lardata/DetectorInfoServices/DetectorClocksService.h" +#include "larcore/Geometry/Geometry.h" +#include "larcore/CoreUtils/ServiceUtil.h" +#include "larcorealg/Geometry/GeometryCore.h" +#include "larsim/MCCheater/BackTrackerService.h" +#include "larsim/MCCheater/ParticleInventoryService.h" +#include "larcore/CoreUtils/ServiceUtil.h" + +class SPAna; + + +class SPAna : public art::EDAnalyzer { +public: + explicit SPAna(fhicl::ParameterSet const& p); + // The compiler-generated destructor is fine for non-base + // classes without bare pointers or other resource use. + + // Plugins should not be copied or assigned. + SPAna(SPAna const&) = delete; + SPAna(SPAna&&) = delete; + SPAna& operator=(SPAna const&) = delete; + SPAna& operator=(SPAna&&) = delete; + + // Required functions. + void analyze(art::Event const& e) override; + + void Clear(); + void setBranches(); + void setMetaBranches(TTree *); + + sbn::HitInfo MakeHit(const recob::Hit &hit, + unsigned hkey, + const art::Ptr &sp, + const geo::GeometryCore *geo, + const detinfo::DetectorClocksData &dclock, + const cheat::BackTrackerService *bt_serv); + + std::pair, std::vector> ParticleTrueHits( + const simb::MCParticle &particle, + const std::map>> id_to_ide_map, + const detinfo::DetectorPropertiesData &dprop, + const geo::GeometryCore *geo); + + void BuildHitRanges(const std::vector> &hits); + void BuildWireRanges(const std::vector> &wires); + + bool FindHit(unsigned channel, unsigned tick); + bool FindWire(unsigned channel, unsigned tick); + + void respondToOpenInputFile(const art::FileBlock& fb) override { + (void) fb; + _fileno ++; + } + +private: + std::vector fHitProducers; + std::vector fWireProducers; + art::InputTag fG4Producer; + art::InputTag fSimChannelProducer; + + TTree *tRecoHits; + TTree *tRecoWire; + TTree *tTrueHits; + + int _fileno; + int _run; + int _subrun; + int _evt; + int _cryo; + + std::vector _reco_hits; + std::vector _true_hits; + + std::vector _true_hit_dir_x; + std::vector _true_hit_dir_y; + std::vector _true_hit_dir_z; + std::vector _true_hits_has_hit; + std::vector _true_hits_has_wire; + + std::vector _reco_wire_start; + std::vector _reco_wire_end; + std::vector _reco_wire_wire; + std::vector _reco_wire_plane; + std::vector _reco_wire_channel; + std::vector _reco_wire_charge; + std::vector _reco_wire_true_charge; + std::vector _reco_wire_true_energy; + + std::map>> fHitRanges; + std::map>> fWireRanges; +}; + + +SPAna::SPAna(fhicl::ParameterSet const& p) + : EDAnalyzer{p} // , + // More initializers here. +{ + // Call appropriate consumes<>() for any products to be retrieved by this module. + art::ServiceHandle tfs; + tRecoHits = tfs->make("RecoHits", "Reconstructed hits"); + tRecoWire = tfs->make("RecoWire", "Reconstructed wires"); + tTrueHits = tfs->make("TrueHits", "True hits"); + + setBranches(); + + _fileno = 0; + _cryo = p.get("Cryostat"); + + fHitProducers = p.get>("HitProducers"); + fWireProducers = p.get>("WireProducers"); + fG4Producer = p.get("G4Producer"); + fSimChannelProducer = p.get("SimChannelProducer"); +} + +void SPAna::BuildHitRanges(const std::vector> &hits) { + fHitRanges.clear(); + for (const art::Ptr &h: hits) fHitRanges[h->Channel()].push_back({h->StartTick(), h->EndTick()}); +} + +void SPAna::BuildWireRanges(const std::vector> &wires) { + fWireRanges.clear(); + for (const art::Ptr &w: wires) { + unsigned channel = w->Channel(); + for (auto const &range: w->SignalROI().get_ranges()) { + fWireRanges[channel].push_back({range.begin_index(), range.end_index()}); + } + } +} + +bool FindRange(const std::map>> &ranges, unsigned channel, unsigned tick) { + if (!ranges.count(channel)) return false; + + const std::vector> &thisrange = ranges.at(channel); + for (const std::pair &r: thisrange) { + if (tick >= r.first && tick <= r.second) return true; + } + + return false; +} + +bool SPAna::FindHit(unsigned channel, unsigned tick) { + return FindRange(fHitRanges, channel, tick); +} + +bool SPAna::FindWire(unsigned channel, unsigned tick) { + return FindRange(fWireRanges, channel, tick); +} + +void SPAna::setBranches() { + setMetaBranches(tRecoHits); + setMetaBranches(tRecoWire); + setMetaBranches(tTrueHits); + + tRecoHits->Branch("h", &_reco_hits); + + tTrueHits->Branch("t", &_true_hits); + tTrueHits->Branch("t.dir_x", &_true_hit_dir_x); + tTrueHits->Branch("t.dir_y", &_true_hit_dir_y); + tTrueHits->Branch("t.dir_z", &_true_hit_dir_z); + tTrueHits->Branch("t.has_hit", &_true_hits_has_hit); + tTrueHits->Branch("t.has_wire", &_true_hits_has_wire); + + tRecoWire->Branch("w.start", &_reco_wire_start); + tRecoWire->Branch("w.end", &_reco_wire_end); + tRecoWire->Branch("w.wire", &_reco_wire_wire); + tRecoWire->Branch("w.plane", &_reco_wire_plane); + tRecoWire->Branch("w.channel", &_reco_wire_channel); + tRecoWire->Branch("w.charge", &_reco_wire_charge); + tRecoWire->Branch("w.true_charge", &_reco_wire_true_charge); + tRecoWire->Branch("w.true_energy", &_reco_wire_true_energy); + +} + +void SPAna::setMetaBranches(TTree *t) { + t->Branch("fileno", &_fileno, "fileno/i"); + t->Branch("run", &_run, "run/i"); + t->Branch("sub", &_subrun, "sub/i"); + t->Branch("evt", &_evt, "evt/i"); + t->Branch("cryo", &_cryo, "cryo/i"); +} + +void SPAna::Clear() { + _reco_hits.clear(); + _true_hits.clear(); + + _true_hit_dir_x.clear(); + _true_hit_dir_y.clear(); + _true_hit_dir_z.clear(); + _true_hits_has_hit.clear(); + _true_hits_has_wire.clear(); + + _reco_wire_start.clear(); + _reco_wire_end.clear(); + _reco_wire_wire.clear(); + _reco_wire_plane.clear(); + _reco_wire_channel.clear(); + _reco_wire_charge.clear(); + _reco_wire_true_charge.clear(); + _reco_wire_true_energy.clear(); +} + +std::map>> PrepSimChannels(const std::vector> &simchannels, const geo::GeometryCore &geo) { + std::map>> ret; + + for (const art::Ptr sc : simchannels) { + // Lookup the wire of this channel + raw::ChannelID_t channel = sc->Channel(); + std::vector maybewire = geo.ChannelToWire(channel); + geo::WireID thisWire; // Default constructor makes invalid wire + if (maybewire.size()) thisWire = maybewire[0]; + + for (const auto &item : sc->TDCIDEMap()) { + for (const sim::IDE &ide: item.second) { + // indexing initializes empty vector + ret[abs(ide.trackID)].push_back({thisWire, item.first, &ide}); + } + } + } + return ret; +} + +// Turn a particle position to a space-charge induced position +geo::Point_t TrajectoryToWirePosition(const geo::Point_t &loc, const geo::TPCID &tpc) { + auto const* sce = lar::providerFrom(); + art::ServiceHandle geom; + + geo::Point_t ret = loc; + + // Returned X is the drift -- multiply by the drift direction to undo this + int corr = geom->TPC(tpc).DriftDir().X(); + + if (sce && sce->EnableSimSpatialSCE()) { + geo::Vector_t offset = sce->GetPosOffsets(ret); + + ret.SetX(ret.X() + corr * offset.X()); + ret.SetY(ret.Y() + offset.Y()); + ret.SetZ(ret.Z() + offset.Z()); + } + + return ret; +} + +// Turn a space-charge induced position to a trajectory Position +geo::Point_t WireToTrajectoryPosition(const geo::Point_t &loc, const geo::TPCID &tpc) { + auto const* sce = lar::providerFrom(); + + geo::Point_t ret = loc; + + if (sce && sce->EnableSimSpatialSCE()) { + geo::Vector_t offset = sce->GetCalPosOffsets(ret, tpc.TPC); + + ret.SetX(ret.X() + offset.X()); + ret.SetY(ret.Y() + offset.Y()); + ret.SetZ(ret.Z() + offset.Z()); + } + + return ret; + +} + +std::pair, std::vector> SPAna::ParticleTrueHits( + const simb::MCParticle &particle, + const std::map>> id_to_ide_map, + const detinfo::DetectorPropertiesData &dprop, + const geo::GeometryCore *geo) { + + // Organize deposition info into per-wire true "Hits" -- key is the Channel Number + std::map truehits; + + const std::vector> empty; + const std::vector> &particle_ides = id_to_ide_map.count(particle.TrackId()) ? id_to_ide_map.at(particle.TrackId()) : empty; + + for (auto const &ide_tup: particle_ides) { + const geo::WireID &w = std::get<0>(ide_tup); + short tick = std::get<1>(ide_tup); + + // cut on cryostat + if ((int)w.Cryostat != _cryo) continue; + + unsigned c = geo->PlaneWireToChannel(w); + const sim::IDE *ide = std::get<2>(ide_tup); + + // Set stuff + truehits[c].cryo = w.Cryostat; + truehits[c].tpc = w.TPC; + truehits[c].plane = w.Plane; + truehits[c].wire = w.Wire; + truehits[c].channel = c; + + // Average stuff using charge-weighting + float old_elec = truehits[c].nelec; + float new_elec = old_elec + ide->numElectrons; + truehits[c].p.x = (truehits[c].p.x*old_elec + ide->x*ide->numElectrons) / new_elec; + truehits[c].p.y = (truehits[c].p.y*old_elec + ide->y*ide->numElectrons) / new_elec; + truehits[c].p.z = (truehits[c].p.z*old_elec + ide->z*ide->numElectrons) / new_elec; + truehits[c].time = (truehits[c].time*old_elec + tick*ide->numElectrons) / new_elec; + + // Also get the position with space charge un-done + geo::Point_t ide_p(ide->x, ide->y, ide->z); + geo::Point_t ide_p_scecorr = WireToTrajectoryPosition(ide_p, w); + + truehits[c].p_scecorr.x = (truehits[c].p_scecorr.x*old_elec + ide_p_scecorr.x()*ide->numElectrons) / new_elec; + truehits[c].p_scecorr.y = (truehits[c].p_scecorr.y*old_elec + ide_p_scecorr.y()*ide->numElectrons) / new_elec; + truehits[c].p_scecorr.z = (truehits[c].p_scecorr.z*old_elec + ide_p_scecorr.z()*ide->numElectrons) / new_elec; + + // Sum stuff + truehits[c].nelec += ide->numElectrons; + truehits[c].e += ide->energy; + truehits[c].ndep += 1; + } + + // Compute widths + for (auto const &ide_tup: particle_ides) { + const geo::WireID &w = std::get<0>(ide_tup); + + // cut on cryostat + if ((int)w.Cryostat != _cryo) continue; + + unsigned c = geo->PlaneWireToChannel(w); + const sim::IDE *ide = std::get<2>(ide_tup); + + geo::Point_t ide_p(ide->x, ide->y, ide->z); + geo::Point_t ide_p_scecorr = WireToTrajectoryPosition(ide_p, w); + + // Average stuff using charge-weighting + float this_elec = ide->numElectrons; + + truehits[c].p_width.x += (ide_p.x() - truehits[c].p.x) * (ide_p.x() - truehits[c].p.x) * this_elec / truehits[c].nelec; + truehits[c].p_width.y += (ide_p.y() - truehits[c].p.y) * (ide_p.y() - truehits[c].p.y) * this_elec / truehits[c].nelec; + truehits[c].p_width.z += (ide_p.z() - truehits[c].p.z) * (ide_p.z() - truehits[c].p.z) * this_elec / truehits[c].nelec; + + truehits[c].p_scecorr_width.x += (ide_p_scecorr.x() - truehits[c].p_scecorr.x) * (ide_p_scecorr.x() - truehits[c].p_scecorr.x) * this_elec / truehits[c].nelec; + truehits[c].p_scecorr_width.y += (ide_p_scecorr.y() - truehits[c].p_scecorr.y) * (ide_p_scecorr.y() - truehits[c].p_scecorr.y) * this_elec / truehits[c].nelec; + truehits[c].p_scecorr_width.z += (ide_p_scecorr.z() - truehits[c].p_scecorr.z) * (ide_p_scecorr.z() - truehits[c].p_scecorr.z) * this_elec / truehits[c].nelec; + } + + // Convert to vector + std::vector truehits_v; + for (auto const &p: truehits) { + truehits_v.push_back(p.second); + } + + // Save true directions + std::vector truehitdirs; + + // Compute the time of each hit + for (sbn::TrueHit &h: truehits_v) { + // TODO: fix magic number + h.time -= 2900; // == (G4RefTime - TriggerOffsetTPC)/TickPeriod = (1500 - 340)/0.4 + double xdrift = abs(h.p.x - geo->Plane(geo::PlaneID(h.cryo, h.tpc, 0)).GetCenter().X()); + h.tdrift = xdrift / dprop.DriftVelocity(); + } + + // Compute the pitch of each hit and order it in the trajectory + for (sbn::TrueHit &h: truehits_v) { + // Use the SCE-undone hit since this matches to the Trajectory + TVector3 h_p(h.p_scecorr.x, h.p_scecorr.y, h.p_scecorr.z); + + TVector3 direction; + float closest_dist = -1.; + int traj_index = -1; + for (unsigned i_traj = 0; i_traj < particle.NumberTrajectoryPoints(); i_traj++) { + if (closest_dist < 0. || (particle.Position(i_traj).Vect() - h_p).Mag() < closest_dist) { + direction = particle.Momentum(i_traj).Vect().Unit(); + closest_dist = (particle.Position(i_traj).Vect() - h_p).Mag(); + traj_index = i_traj; + } + } + + sbn::Vector3D dir_v; + dir_v.x = direction.x(); + dir_v.y = direction.y(); + dir_v.z = direction.z(); + truehitdirs.push_back(dir_v); + + // If we got a direction, get the pitch + if (closest_dist >= 0. && direction.Mag() > 1e-4) { + geo::PlaneID plane(h.cryo, h.tpc, h.plane); + float angletovert = geo->WireAngleToVertical(geo->View(plane), plane) - 0.5*::util::pi<>(); + float cosgamma = abs(cos(angletovert) * direction.Z() + sin(angletovert) * direction.Y()); + float pitch = geo->WirePitch(plane) / cosgamma; + h.pitch = pitch; + } + else { + h.pitch = -1.; + } + // And the pitch induced by SCE + if (closest_dist >= 0. && direction.Mag() > 1e-4) { + geo::PlaneID plane(h.cryo, h.tpc, h.plane); + float angletovert = geo->WireAngleToVertical(geo->View(plane), plane) - 0.5*::util::pi<>(); + + TVector3 loc_mdx_v = h_p - direction * (geo->WirePitch(geo->View(plane)) / 2.); + TVector3 loc_pdx_v = h_p + direction * (geo->WirePitch(geo->View(plane)) / 2.); + + // Convert types for helper functions + geo::Point_t loc_mdx(loc_mdx_v.X(), loc_mdx_v.Y(), loc_mdx_v.Z()); + geo::Point_t loc_pdx(loc_pdx_v.X(), loc_pdx_v.Y(), loc_pdx_v.Z()); + geo::Point_t h_p_point(h_p.X(), h_p.Y(), h_p.Z()); + + loc_mdx = TrajectoryToWirePosition(loc_mdx, plane); + loc_pdx = TrajectoryToWirePosition(loc_pdx, plane); + + // Direction at wires + geo::Vector_t dir = (loc_pdx - loc_mdx) / (loc_mdx - loc_pdx).r(); + + // Pitch at wires + double cosgamma = std::abs(std::sin(angletovert)*dir.Y() + std::cos(angletovert)*dir.Z()); + double pitch; + if (cosgamma) { + pitch = geo->WirePitch(geo->View(plane))/cosgamma; + } + else { + pitch = 0.; + } + + // Now bring that back to the particle trajectory + geo::Point_t loc_w = TrajectoryToWirePosition(h_p_point, plane); + + geo::Point_t locw_pdx_traj = WireToTrajectoryPosition(loc_w + pitch*dir, plane); + geo::Point_t loc = WireToTrajectoryPosition(loc_w, plane); + + h.pitch_sce = (locw_pdx_traj - loc).R(); + } + else { + h.pitch_sce = -1.; + } + + // And the trajectory location + h.itraj = traj_index; + + // And the residual range of the hit + h.rr = 0.; + if (traj_index >= 0) { + for (int i_traj = traj_index+1; i_traj < (int)particle.NumberTrajectoryPoints(); i_traj++) { + h.rr += (particle.Position(i_traj).Vect() - particle.Position(i_traj-1).Vect()).Mag(); + } + + // Also account for the distance from the Hit point to the matched trajectory point + double hit_distance_along_particle = (h_p - particle.Position(traj_index).Vect()).Dot(particle.Momentum(traj_index).Vect().Unit()); + h.rr += -hit_distance_along_particle; + } + } + + return {truehits_v, truehitdirs}; + +} + +sbn::HitInfo SPAna::MakeHit(const recob::Hit &hit, + unsigned hkey, + const art::Ptr &sp, + const geo::GeometryCore *geo, + const detinfo::DetectorClocksData &dclock, + const cheat::BackTrackerService *bt_serv) { + + // TrackHitInfo to save + sbn::HitInfo h; + + // information from the hit object + h.integral = hit.Integral(); + h.sumadc = hit.ROISummedADC(); + h.width = hit.RMS(); + h.time = hit.PeakTime(); + h.mult = hit.Multiplicity(); + h.wire = hit.WireID().Wire; + h.plane = hit.WireID().Plane; + h.channel = geo->PlaneWireToChannel(hit.WireID()); + h.tpc = hit.WireID().TPC; + h.end = hit.EndTick(); + h.start = hit.StartTick(); + h.id = (int)hkey; + + // Do back-tracking on each hit + if (bt_serv) { + // The default BackTracking function goes from (peak - width, peak + width). + // + // This time range does not match well hits with a non-Gaussian shape where + // the Gaussian-fit-width does not replicate the width of the pulse. + // + // Instead, we use the Hit (start, end) time range. This is also more relevant + // for (e.g.) the SummedADC charge extraction method. + // + // Don't use this: + // std::vector ides = bt_serv->HitToTrackIDEs(dclock, hit); + // + // Use this: + std::vector ides = bt_serv->ChannelToTrackIDEs(dclock, hit.Channel(), hit.StartTick(), hit.EndTick()); + + h.truth.e = 0.; + h.truth.nelec = 0.; + + for (const sim::TrackIDE &ide: ides) { + h.truth.e += ide.energy; + h.truth.nelec += ide.numElectrons; + } + } + else { + h.truth.e = -1.; + h.truth.nelec = -1.; + } + + + // Save SpacePoint information + if (sp) { + h.sp.x = sp->position().x(); + h.sp.y = sp->position().y(); + h.sp.z = sp->position().z(); + + h.hasSP = true; + } + else { + h.hasSP = false; + } + + return h; +} + +void SPAna::analyze(art::Event const& e) +{ + unsigned evt = e.event(); + unsigned sub = e.subRun(); + unsigned run = e.run(); + _evt = evt; + _subrun = sub; + _run = run; + + std::cout << "[SPAna::analyze] Run: " << run << ", SubRun: " << sub << ", Event: "<< evt << ", Is Data: " << e.isRealData() << std::endl; + + Clear(); + + // Load services + const geo::GeometryCore *geometry = lar::providerFrom(); + auto const clock_data = art::ServiceHandle()->DataFor(e); + auto const dprop = + art::ServiceHandle()->DataFor(e, clock_data); + + art::ServiceHandle bt_serv; + + + // Collect products + + // Hits + std::vector> hitList; + for (const art::InputTag &t: fHitProducers) { + art::ValidHandle> hitHandle = e.getValidHandle>(t); + art::fill_ptr_vector(hitList, hitHandle); + } + // Add spacepoints later? + // art::FindManyP allHitSPs(hitList, e, fPFPproducer); + + // Wires + std::vector> wireList; + for (const art::InputTag &t: fWireProducers) { + art::ValidHandle> wireHandle = e.getValidHandle>(t); + art::fill_ptr_vector(wireList, wireHandle); + } + + // Truth + std::vector> mcparticles; + art::ValidHandle> mcparticle_handle = e.getValidHandle>(fG4Producer); + art::fill_ptr_vector(mcparticles, mcparticle_handle); + + std::vector> simchannels; + art::ValidHandle> simchannel_handle = e.getValidHandle>(fSimChannelProducer); + art::fill_ptr_vector(simchannels, simchannel_handle); + + // Prep matching info + std::map>> id_to_ide_map = PrepSimChannels(simchannels, *geometry); + BuildHitRanges(hitList); + BuildWireRanges(wireList); + + // Save hits + for (unsigned i = 0; i < hitList.size(); i++) { + _reco_hits.push_back(MakeHit(*hitList[i], hitList[i].key(), {} /*allHitSps.at(hitList[i].key())*/, geometry, clock_data, bt_serv.get())); + } + + // Save wires + for (unsigned i = 0; i < wireList.size(); i++) { + art::Ptr wire = wireList[i]; + + unsigned channel = wire->Channel(); + unsigned plane_id = geometry->ChannelToWire(wire->Channel()).at(0).Plane; + unsigned wire_id = geometry->ChannelToWire(wire->Channel()).at(0).Wire; + for (auto const &range: wire->SignalROI().get_ranges()) { + + _reco_wire_start.push_back(range.begin_index()); + _reco_wire_end.push_back(range.end_index()); + _reco_wire_wire.push_back(wire_id); + _reco_wire_plane.push_back(plane_id); + _reco_wire_channel.push_back(channel); + + float charge = 0; + for (float val: range) charge += val; + _reco_wire_charge.push_back(charge); + + float true_charge = -1; + float true_energy = -1; + // Do back-tracking on each wire + if (bt_serv.get()) { + std::vector ides = bt_serv->ChannelToTrackIDEs(clock_data, channel, range.begin_index(), range.end_index()); + + true_charge = 0; + true_energy = 0; + + for (const sim::TrackIDE &ide: ides) { + true_energy += ide.energy; + true_charge += ide.numElectrons; + } + } + _reco_wire_true_charge.push_back(true_charge); + _reco_wire_true_energy.push_back(true_energy); + } + } + + // Save truth info + for (unsigned i = 0; i < mcparticles.size(); i++) { + const simb::MCParticle &p = *mcparticles[i]; + auto ret = ParticleTrueHits(p, id_to_ide_map, dprop, geometry); + + for (const sbn::TrueHit &th: ret.first) _true_hits.push_back(th); + for (const sbn::Vector3D &d: ret.second) { + _true_hit_dir_x.push_back(d.x); + _true_hit_dir_y.push_back(d.y); + _true_hit_dir_z.push_back(d.z); + } + } + + // true to reco matching + for (const sbn::TrueHit &th: _true_hits) { + _true_hits_has_hit.push_back((int)FindHit(th.channel, th.time)); + _true_hits_has_wire.push_back((int)FindWire(th.channel, th.time)); + } + + tRecoHits->Fill(); + tTrueHits->Fill(); + tRecoWire->Fill(); + +} + +DEFINE_ART_MODULE(SPAna) diff --git a/icaruscode/TPC/SPAna/run_spana_icarus.fcl b/icaruscode/TPC/SPAna/run_spana_icarus.fcl new file mode 100644 index 000000000..9224fe419 --- /dev/null +++ b/icaruscode/TPC/SPAna/run_spana_icarus.fcl @@ -0,0 +1,67 @@ +#include "services_common_icarus.fcl" +#include "simulationservices_icarus.fcl" + +#include "larproperties.fcl" +#include "backtrackerservice.fcl" +#include "particleinventoryservice.fcl" +#include "detectorproperties_icarus.fcl" +#include "spacecharge.fcl" + +#include "spana_icarus.fcl" + +process_name: SPAna + +services: +{ + @table::icarus_basic_services + @table::icarus_wirecalibration_services + @table::icarus_backtracking_services # from `simulationservices_icarus.fcl` + SpaceChargeService: @local::icarus_spacecharge +} +services.SpaceChargeService: { + EnableCalEfieldSCE: false + EnableCalSpatialSCE: false + EnableCorrSCE: false + EnableSimEfieldSCE: false + EnableSimSpatialSCE: false + InputFilename: "SCEoffsets/SCEoffsets_ICARUS_E500_voxelTH3.root" + RepresentationType: "Voxelized_TH3" + service_provider: "SpaceChargeServiceICARUS" +} + +services.TFileService.fileName: "SPAna.root" + +#Start each new event with an empty event. +source: +{ + module_type: RootInput + maxEvents: -1 # Number of events to create +} + +outputs: +{ +} + + +# Define and configure some modules to do work on each event. +# First modules are defined; they are scheduled later. +# Modules are grouped by type. +physics: +{ + producers:{} + + filters: {} + + # evtfilter: [filter] + + analyzers: + { + spanaE: @local::spana2d_east + spanaW: @local::spana2d_west + } + + runprod: [] + stream1: [spanaE, spanaW] + trigger_paths: [] + end_paths: [ stream1 ] +} diff --git a/icaruscode/TPC/SPAna/spana_icarus.fcl b/icaruscode/TPC/SPAna/spana_icarus.fcl new file mode 100644 index 000000000..cebb35042 --- /dev/null +++ b/icaruscode/TPC/SPAna/spana_icarus.fcl @@ -0,0 +1,39 @@ +BEGIN_PROLOG + +spana1d_east: { + module_type: SPAna + Cryostat: 0 + HitProducers: ["gaushit1dTPCEE", "gaushit1dTPCEW"] + WireProducers: ["roifinder1d:PHYSCRATEDATATPCEE", "roifinder1d:PHYSCRATEDATATPCEW"] + G4Producer: "largeant" + SimChannelProducer: "daq:simpleSC" +} + +spana1d_west: { + module_type: SPAna + Cryostat: 1 + HitProducers: ["gaushit1dTPCWE", "gaushit1dTPCWW"] + WireProducers: ["roifinder1d:PHYSCRATEDATATPCWE", "roifinder1d:PHYSCRATEDATATPCWW"] + G4Producer: "largeant" + SimChannelProducer: "daq:simpleSC" +} + +spana2d_east: { + module_type: SPAna + Cryostat: 0 + HitProducers: ["gaushit2dTPCEE", "gaushit2dTPCEW"] + WireProducers: ["decon2droiEE:gauss", "decon2droiEW:gauss"] + G4Producer: "largeant" + SimChannelProducer: "daq:simpleSC" +} + +spana2d_west: { + module_type: SPAna + Cryostat: 1 + HitProducers: ["gaushit2dTPCWE", "gaushit2dTPCWW"] + WireProducers: ["decon2droiWE:gauss", "decon2droiWW:gauss"] + G4Producer: "largeant" + SimChannelProducer: "daq:simpleSC" +} + +END_PROLOG diff --git a/icaruscode/TPC/SignalProcessing/HitFinder/hitfindermodules_icarus.fcl b/icaruscode/TPC/SignalProcessing/HitFinder/hitfindermodules_icarus.fcl index 722f6355c..bda1f1dd9 100644 --- a/icaruscode/TPC/SignalProcessing/HitFinder/hitfindermodules_icarus.fcl +++ b/icaruscode/TPC/SignalProcessing/HitFinder/hitfindermodules_icarus.fcl @@ -100,7 +100,6 @@ gaus_hitfinder_icarus.HitFinderToolVec.CandidateHitsPlane2: gaus_hitfinder_icarus.HitFinderToolVec.CandidateHitsPlane2.Plane: 2 gaus_hitfinder_icarus.HitFinderToolVec.CandidateHitsPlane2.RoiThreshold: 9. - # Define icarus version of test version gausshit finder gauss_hitfinder_icarus: @local::gaus_hitfinder