aboutsummaryrefslogtreecommitdiff
path: root/utils
AgeCommit message (Collapse)Author
2021-01-05print latex tabletlatorre
2021-01-05update dc flasher and noise codetlatorre
2021-01-05add --fast argument to dm-searchtlatorre
2021-01-05add --fasttlatorre
2021-01-05hack to get rid of flasher and muon events in breakdown sampletlatorre
2021-01-04update plot-dctlatorre
2021-01-04get rid of nhit_threshtlatorre
2021-01-04update radius cut in dctlatorre
2021-01-03add --universe argument to dm-searchtlatorre
2021-01-03speed up submit-grid-jobstlatorre
This commit updates submit-grid-jobs to not run zdab-cat if the file has already been submitted.
2021-01-03add a numba optimized version of interptlatorre
2021-01-03speed up read_mcpl()tlatorre
2021-01-03reduce memory usage by creating weights dict earlytlatorre
2021-01-03cache results from get_events()tlatorre
2020-12-27reduce memory usagetlatorre
2020-12-27reduce memory usagetlatorre
2020-12-27reduce memory usagetlatorre
2020-12-27add dm-search to Makefiletlatorre
2020-12-25don't exit if some runs don't have MCtlatorre
2020-12-25merge lefttlatorre
2020-12-25don't warn about rhdr banks and 50 MHz clock jumps for MCtlatorre
2020-12-22move check earliertlatorre
2020-12-22fix typo from last committlatorre
2020-12-22go back to checking rhdr instead of evtlatorre
2020-12-21fix another rhdr checktlatorre
2020-12-21fix rhdr checktlatorre
2020-12-21only concat files with at least 1 eventtlatorre
2020-12-21use correct bins when calculation thresholdtlatorre
2020-12-21update dm-search to only look for mediator masses up to 1 GeVtlatorre
2020-12-21update bins for muon histogramstlatorre
This commit updates the bins for the muon histograms to not go too far above 2 GeV. The reason is that at these energies most muons will exit the detetor and my energy reconstruction doesn't work too well. I also updated chi2 and dm-search to use the flux_weight when sampling from the MC.
2020-12-16use a hash to merge weights with MC datatlatorre
2020-12-16print out livetime from run_info.log in dm-searchtlatorre
2020-12-16python 3 fixtlatorre
2020-12-15add code to reweight the tau neutrino eventstlatorre
This commit updates the code to reweight the MC data from tau neutrinos since I stupidly simulated the muon neutrino flux instead of the tau neutrino flux.
2020-12-09update discovery threshold againtlatorre
2020-12-09ensure discovery threshold is at least 0tlatorre
2020-12-09speed up chi2 and dm-searchtlatorre
2020-12-09use nanpercentile in get_events()tlatorre
2020-12-09discard negative weightstlatorre
2020-12-09fix bug introduced in previous committlatorre
This commit fixes a bug I introduced earlier in chi2 and dm-search since we want to remove runs not in the MC for both the signal and atmospheric samples.
2020-12-09fix bug when running chi2 --pulltlatorre
2020-12-09update logger.py to be python3 compatibletlatorre
This commite updates the Logger class to use str instead of basestring.
2020-12-09don't delete mc data from runs not in data when running chi2 --coveragetlatorre
2020-12-08speed up dm-searchtlatorre
2020-12-08use weighted MC when running dm-search --testtlatorre
2020-12-08fix more typostlatorre
2020-12-08fix another typotlatorre
2020-12-08fix typo from last committlatorre
2020-12-08update dm-search and chi2 to refit with GENIE systematicstlatorre
2020-12-08don't sample DM events when running --testtlatorre
his program. If not, see <https://www.gnu.org/licenses/>. """ Script to reprocess ZDAB files from SNO and write out orphans. Copied mostly from zdab-reprocess, but I changed the command file template. To reprocess a single file: $ ./zdab-reprocess-orphans FILENAME and to batch reprocess: $ ./zdab-reprocessi-orphans FILENAME FILENAME ... By default, the reprocessed files will be stored in the current working directory with _reprocessed appended to the name, but this can be changed by passing a different suffix on the command line. For example: $ ./zdab-reprocess-orphans --suffix [suffix] SNOCR_00000100004_000_p2.xzdab will produce a file named SNOCR_00000100004_000_p2_suffix.xzdab. """ from __future__ import print_function, division import string import os import re def GetSettingsFile(run): """ Returns the settings file for a given run number. Taken from GetSettingsFile() in autosno/lib/AutoSNO_Module.pm. """ if run < 19999: return "load_d2o_settings.cmd" elif run < 33907: return "load_salt_settings.cmd" elif run < 39999: return "load_d2o_2_settings.cmd" elif run < 67900: return "load_ncd_settings.cmd" else: return "load_h2o_settings.cmd" # from autosno/run_lists/RunList_clock_force_no_fix.dat RUN_LIST_CLOCK_FORCE_NO_FIX = [ 11371, 11377, 11399, 11436, 11575, 11681, 11701, 11911, 11976, 12187, 12233, 13401, 13423, 13431, 13874, 14264, 51456, ] def get_clock_force_no_fix_string(run): if run in RUN_LIST_CLOCK_FORCE_NO_FIX: return "set bank TUPK 1 word 3 to -1" else: return "" # This template comes from autosno/prod/Reconstruct.pm. I also added the lines # which autosno adds in Run() of the Reconstruct module. SNOMAN_TEMPLATE = \ """ * This file was automatically created by zdab-reprocess. * * Note: This command file was created by attempting to "reverse engineer" how * autosno fits neutrino runs. * * Author: Anthony LaTorre * Created on August 26, 2019. file inp 1 @input file out 1 @output titles @dqxx titles @anxx_pca titles @anxx_neutrino @@@settings_file $flt_set_prescale @prescale_setting $define_test 1 $line_1 'nhits float_equals EV+$KEV_NPM:@lower_nhit @upper_nhit;' $starting_seed @seed $starting_seed_2 @seed2 @clock_force_no_fix * We want to process the entire file. $num_events 0 * Titles and database checks. set bank TMTT 1 word 0 to 1 * HCA Control. $hca_mode $perform_correction $hca_correction $new_correction $hca_consecutive $yes * QRC Control. * We need the cross talk cut. $apply_xtalk_cut *Settings for muon fitters. * Set the charge threshold for multi photon timing correction $set_mpca_threshold 1.5 titles mpca_data.dat * Turn on the channel dependent gain calibration set bank TCAL 1 word 4 to 54 *Load charge probability titles file titles charge_photoelectrons.dat * Low nhit prescale $flt_set_prescale_2 0.0 * Primary NHIT Cut. $enable_test 1 $define_test 1 $line_1 'nhits float_equals EV+$KEV_NPM:18.5 9999.;' * Primary Fitter NHIT Cut. $enable_test 3 $define_test 3 $line_1 'nhits float_equals EV+$KEV_NPM:9.5 999.;' * NCD Event Tag. $enable_test 4 $define_test 4 $line_1 'ncd float_equals EV+$KEV_NCD_STATUS:0.5 3.5;' * Muon Candidate Tag. $enable_test 5 $define_test 5 $line_1 'nhits float_equals EV+$KEV_NPM:249.5 9999.;' * A do nothing filter that allows autosno to filter maintenece events. $enable_test 9 $define_test 9 $line_1 'one equals 1.0;' define event_loop call INP call UPK call CAL call FLT($PERM_BANK_CUT) if_not_ok goto WRITE_OUT EVENT_REJECT: call FLT(9) * User defined cut. Default passes all. if_not_ok quit_event call FLT($JUNK_CUT) * Junk. We still want to keep NCD events. if_not_ok goto WRITE_OUT quit_event WRITE_OUT: call PRU call OUT quit_event end_def *Pruning Control $prune_pix $drop $prune_px $drop $prune_pf $drop $prune_nesg $drop $prune_nes $drop $prune_nemg $drop * Run with the Database @@run_snodb !endfile member=fit_neutrino """.strip() class MyTemplate(string.Template): delimiter = '@' def splitext(path): """ Like os.path.splitext() except it returns the full extension if the filename has multiple extensions, for example: splitext('foo.tar.gz') -> 'foo', '.tar.gz' """ full_root, full_ext = os.path.splitext(path) while True: root, ext = os.path.splitext(full_root) if ext: full_ext = ext + full_ext full_root = root else: break return full_root, full_ext def GetRunRange(run, num=1000): """ Returns the run range string for ANXX neutrino files. Example: GetRunRange(10000,1000) -> "10000-10999" Comes from GetRunRange() in autosno/lib/FileUtil.pm. """ base_run = int(run//num)*num return "%i-%i" % (base_run, base_run + num - 1) def get_anxx_neutrino_filename(run, dir): """ Returns the ANXX neutrino file with the highest pass number for a given run. The files should all be located in dir. FIXME: Is it OK to just get the file with the highest pass number? Technically this info should come from autosno, but it's not clear exactly how it's getting the pass number in the function GetPassNum() in AutoSNO_Module.pm. """ filename_re = "anxx_nu_(\d+)_p(\d+).dat" p = re.compile(filename_re) filenames = [] for filename in os.listdir(dir): match = p.match(filename) if int(match.group(1)) != run: continue filenames.append((match.group(2),filename)) return sorted(filenames)[-1][1] def get_anxx_pca_next(run,dir): """ Returns the next ANXX PCA file after `run`. See Get_ANxx_PCA_Next() in autosno/lib/GetTitles.pm. """ filename_re = "anxx_pca_(\d+)_p(\d+).dat" p = re.compile(filename_re) filenames = [] for filename in os.listdir(dir): match = p.match(filename) if match is None: continue if int(match.group(1)) < run: continue filenames.append((match.group(1),match.group(2),filename)) return sorted(filenames)[-1][2] if __name__ == '__main__': import argparse import sys import os from subprocess import check_call parser = argparse.ArgumentParser("reprocess zdabs") parser.add_argument("filenames", nargs="+", help="filenames of zdabs") parser.add_argument("--dir", default=None, help="directory to store reprocessed zdab files") parser.add_argument("--lower-nhit", type=float, default=18.5, help="lower value for primary nhit cut") parser.add_argument("--upper-nhit", type=float, default=9999, help="upper value for primary nhit cut") parser.add_argument("--suffix", default="orphans", help="suffix appended to output filename") args = parser.parse_args() filename_re = "(?:SNOCR_|SNO)(\d+)_\w+.x?zdab" p = re.compile(filename_re) template = MyTemplate(SNOMAN_TEMPLATE) for filename in args.filenames: head, tail = os.path.split(filename) match = p.match(tail) if not match: print("Unable to parse filename '%s'" % tail) sys.exit(1) try: run = int(match.group(1)) except Exception as e: print("Unable to convert run number to int: %s" % str(e)) sys.exit(1) root, ext = splitext(tail) output = "%s_%s%s" % (root,args.suffix,ext) if args.dir: output = os.path.join(args.dir,output) cmd_filename = "%s.cmd" % root dqxx_filename = "/sno/mcprod/dqxx/DQXX_%010i.dat" % run anxx_neutrino_directory = "/sno/mcprod/anxx/titles/neutrino/%s" % GetRunRange(run) anxx_neutrino_filename = get_anxx_neutrino_filename(run, anxx_neutrino_directory) anxx_neutrino_filename = os.path.join(anxx_neutrino_directory,anxx_neutrino_filename) anxx_pca_directory = "/sno/mcprod/anxx/titles/pca" anxx_pca_filename = get_anxx_pca_next(run, anxx_pca_directory) anxx_pca_filename = os.path.join(anxx_pca_directory,anxx_pca_filename) with open(cmd_filename, "w") as f: f.write(template.safe_substitute(input=filename, output=output, dqxx=dqxx_filename, anxx_neutrino=anxx_neutrino_filename, anxx_pca=anxx_pca_filename, settings_file=GetSettingsFile(run), prescale_setting=1.0, # no idea what this is lower_nhit = args.lower_nhit, upper_nhit = args.upper_nhit, seed=0, seed2=0, clock_force_no_fix=get_clock_force_no_fix_string(run))) check_call(["snoman.exe","-c",cmd_filename])