Skip to content

Commit

Permalink
Merge pull request #19 from nstelter-slac/bringover_all_scripts
Browse files Browse the repository at this point in the history 
Bringover all scritps from the /rix directory, revert exists scripts to orignal versions for now
  • Loading branch information
@nstelter-slac
nstelter-slac authored Feb 19, 2024
2 parents 20dbded + 3f703a2 commit 118bfea
Show file tree
Hide file tree
Showing 98 changed files with 29,521 additions and 987 deletions.
6 changes: 3 additions & 3 deletions .github/workflows/run-tests.yml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,6 @@ jobs:
pip install scipy
pip install statsmodels
- name: Run tests
run: |
pytest tests/test_FitFunctions.py
# - name: Run tests
# run: |
# pytest tests/test_FitFunctions.py
Binary file added calibrationSuite/OffXavierV4.npy
View file
Binary file not shown.
0 scripts/OffXavierV4_2.npy → calibrationSuite/OffXavierV4_2.npy
100755 → 100644
View file
File renamed without changes.
0 scripts/__init__.py → calibrationSuite/__init__.py
View file
File renamed without changes.
40 changes: 21 additions & 19 deletions scripts/ancillaryMethods.py → calibrationSuite/ancillaryMethods.py
View file
Original file line number Diff line number Diff line change
@@ -1,49 +1,51 @@
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import binned_statistic


def makeProfile(x, y, bins, range=None, spread=False):
## NaN for empty bins are suppressed
## using mean root(N) for non-empty bins to calculate 0 var weights
##
## spread=True to return standard deviation instead of standard error

meansObj = binned_statistic(x, [y, y**2], bins=bins, range=range, statistic="mean")
meansObj = binned_statistic(x, [y, y**2], bins=bins, range=range, statistic='mean')
means, means2 = meansObj.statistic
countsObj = binned_statistic(x, [y, y**2], bins=bins, range=(0, 1), statistic="count")
countsObj = binned_statistic(x, y, bins=bins, range=range, statistic='count')
stdObj = binned_statistic(x, y, bins=bins, range=range, statistic='std')
bin_N = countsObj.statistic
yErr = np.sqrt(means2 - means**2)
usefulBins = np.bitwise_and(bin_N>0, ~np.isnan(means))
if bin_N.sum()==0:
##no data
print("no data in profile")
return None, None, None

##yErr = np.sqrt(means2 - means**2)
yErr = stdObj.statistic
if not spread:
root_N = np.sqrt(bin_N)
root_N[root_N == 0] = root_N[root_N > 0].mean()
yErr = yErr / root_N
root_N[root_N==0] = root_N[usefulBins].mean()
yErr = yErr/root_N
##yErr = yErr.clip(0, 6666666.)
# where is means_result defined?? (ignore from ruff for now)
bin_edges = means_result.bin_edges # noqa: F821
bin_centers = (bin_edges[:-1] + bin_edges[1:]) / 2.0
usefulBins = bin_N > 0
return bin_centers[usefulBins], means[usefulBins], yErr[usefulBins]

bin_edges = meansObj.bin_edges
bin_centers = (bin_edges[:-1] + bin_edges[1:])/2.
return bin_centers[usefulBins], means[usefulBins], yErr[usefulBins]

def plotProfile(x, y, yErr):
plt.errorbar(x=x, y=y, yerr=yErr, linestyle="none", marker=".")

plt.errorbar(x=x, y=y, yerr=yErr, linestyle='none', marker='.')

def selectedClusters(clusters, row, col, lowEnerygCut, highEnergyCut, nPixelCut=4, isSquare=1):
pass


def goodClusters(clusters, row, col, nPixelCut=4, isSquare=None):
##print(clusters)
pixelRowCol = np.bitwise_and((clusters[:, :, 1] == row), (clusters[:, :, 2] == col))
pixelRowCol = np.bitwise_and((clusters[:,:,1] == row),
(clusters[:,:,2] == col))
if isSquare is None:
small = clusters[:, :, 3] < nPixelCut
small = clusters[:,:,3]<nPixelCut
else:
small = np.bitwise_and((clusters[:, :, 3] < nPixelCut), (clusters[:, :, 4] == isSquare))
small = np.bitwise_and((clusters[:,:,3]<nPixelCut), (clusters[:,:,4]==isSquare))
return clusters[np.bitwise_and(small, pixelRowCol)]


def getClusterEnergies(clusters):
##print(clusters)
return clusters[:, 0]
224 changes: 105 additions & 119 deletions scripts/basicSuiteScript.py → calibrationSuite/basicSuiteScript.py
View file
Original file line number Diff line number Diff line change
@@ -1,127 +1,121 @@
import argparse
import numpy as np
import os
from rixSuiteConfig import experimentHash

import matplotlib.pyplot as plt
from matplotlib.ticker import AutoMinorLocator
import sys
from calibrationSuite.rixSuiteConfig import *
##from mfxRixSuiteConfig import *

if os.getenv("foo") == "1":
print("psana1")
from psana1Base import PsanaBase
else:
print("psana2")
from psana2Base import PsanaBase

import h5py
from scipy.optimize import curve_fit ## here?
import calibrationSuite.fitFunctions
import calibrationSuite.ancillaryMethods

def sortArrayByList(a, data):
return [x for _, x in sorted(zip(a, data), key=lambda pair: pair[0])]
return [x for _,x in sorted(zip(a, data), key=lambda pair:pair[0])]


class BasicSuiteScript(PsanaBase):
def __init__(self, analysisType="scan"):
def __init__(self, analysisType='scan'):
super().__init__()
##print("in BasicSuiteScript, inheriting from PsanaBase, type is psana%d" %(self.psanaType))

self.gainModes = {"FH": 0, "FM": 1, "FL": 2, "AHL-H": 3, "AML-M": 4, "AHL-L": 5, "AML-L": 6}
self.ePix10k_cameraTypes = {1: "Epix10ka", 4: "Epix10kaQuad", 16: "Epix10ka2M"}
self.gainModes = {"FH":0, "FM":1, "FL":2, "AHL-H":3, "AML-M":4, "AHL-L":5, "AML-L":6}
self.ePix10k_cameraTypes = {1:"Epix10ka", 4:"Epix10kaQuad", 16:"Epix10ka2M"}
self.camera = 0
##self.outputDir = '/sdf/data/lcls/ds/rix/rixx1003721/results/%s/' %(analysisType)
self.outputDir = "../%s/" % (analysisType)
self.outputDir = '../%s/' %(analysisType)
##self.outputDir = '/tmp'

self.className = self.__class__.__name__

try:
self.location = experimentHash["location"]
except Exception:
self.location = experimentHash['location']
except:
pass
try:
self.exp = experimentHash["exp"]
except Exception:
self.exp = experimentHash['exp']
except:
pass
try:
##if True:
self.ROIfileNames = experimentHash["ROIs"]
##if True:
self.ROIfileNames = experimentHash['ROIs']
self.ROIs = []
for f in self.ROIfileNames:
self.ROIs.append(np.load(f + ".npy"))
self.ROIs.append(np.load(f+'.npy'))
try: ## dumb code for compatibility or expectation
self.ROI = self.ROIs[0]
except Exception:
except:
pass
##if False:
except Exception:
except:
print("had trouble finding", self.ROIfileNames)
self.ROI = None
self.ROIs = None
try:
self.singlePixels = experimentHash["singlePixels"]
except Exception:
self.singlePixels = experimentHash['singlePixels']
except:
self.singlePixels = None
try:
self.regionSlice = experimentHash["regionSlice"]
except Exception:
self.regionSlice = experimentHash['regionSlice']
except:
self.regionSlice = None
if self.regionSlice is not None:
self.sliceCoordinates = [
[self.regionSlice[0].start, self.regionSlice[0].stop],
[self.regionSlice[1].start, self.regionSlice[1].stop],
]
self.sliceCoordinates = [[self.regionSlice[0].start,
self.regionSlice[0].stop],
[self.regionSlice[1].start,
self.regionSlice[1].stop],
]
sc = self.sliceCoordinates
self.sliceEdges = [sc[0][1] - sc[0][0], sc[1][1] - sc[1][0]]

self.sliceEdges = [sc[0][1]-sc[0][0], sc[1][1]-sc[1][0]]
try:
self.fluxSource = experimentHash["fluxSource"]
self.fluxSource = experimentHash['fluxSource']
try:
self.fluxChannels = experimentHash["fluxChannels"]
except Exception:
self.fluxChannels = range(8, 16) ## wave8
self.fluxChannels = experimentHash['fluxChannels']
except:
self.fluxChannels = range(8,16) ## wave8
try:
self.fluxSign = experimentHash["fluxSign"]
except Exception:
self.fluxSign = experimentHash['fluxSign']
except:
self.fluxSign = 1
except Exception:
except:
self.fluxSource = None

## for non-120 Hz running
self.nRunCodeEvents = 0
self.nDaqCodeEvents = 0
self.nBeamCodeEvents = 0
self.runCode = 280
self.daqCode = 281
self.beamCode = 283 ## per Matt
##self.beamCode = 281 ## don't see 283...
self.fakeBeamCode = False

parser = argparse.ArgumentParser(
description="Configures calibration suite, overriding experimentHash",
formatter_class=argparse.ArgumentDefaultsHelpFormatter,
)
parser.add_argument("-e", "--exp", help="experiment")
parser.add_argument("-l", "--location", help="hutch location, e.g. MfxEndstation or DetLab")
parser.add_argument("-r", "--run", type=int, help="run")
parser.add_argument("-R", "--runRange", help="run range, format ...")
parser.add_argument("--fivePedestalRun", type=int, help="5 pedestal run")
parser.add_argument("--fakePedestal", type=str, help="fake pedestal file")
parser.add_argument("-c", "--camera", type=int, help="camera.n")
parser.add_argument("-p", "--path", type=str, help="the base path to the output directory")
parser.add_argument("-n", "--nModules", type=int, help="nModules")
parser.add_argument(
"-d", "--detType", type=str, default="", help="Epix100, Epix10ka, Epix10kaQuad, Epix10ka2M, ..."
)
parser.add_argument("--maxNevents", type=int, default="666666", help="max number of events to analyze")
parser.add_argument(
"--skipNevents", type=int, default=0, help="max number of events to skip at the start of each step"
)
parser.add_argument(
"--configScript",
type=str,
default="experimentSuiteConfig.py",
help="name of python config file to load if any",
)
parser.add_argument("--detObj", help='"raw", "calib", "image"')
parser.add_argument("-f", "--file", type=str, help="run analysis only on file")
parser.add_argument("-L", "--label", type=str, help="analysis label")
parser.add_argument("-t", "--threshold", help="threshold (ADU or keV or wave8) depending on --detObj")
parser.add_argument("--fluxCut", type=float, help="minimum flux to be included in analysis")
parser.add_argument(
"--special",
type=str,
help="comma-separated list of special behaviors - maybe this is too lazy.\
E.g. positiveParity,doKazEvents,...",
description='Configures calibration suite, overriding experimentHash',
formatter_class=argparse.ArgumentDefaultsHelpFormatter
)
parser.add_argument('-e', '--exp', help='experiment')
parser.add_argument('-l', '--location', help='hutch location, e.g. MfxEndstation or DetLab')
parser.add_argument('-r', '--run', type=int, help='run')
parser.add_argument('-R', '--runRange', help='run range, format ...')
parser.add_argument('--fivePedestalRun', type=int, help='5 pedestal run')
parser.add_argument('--fakePedestal', type=str, help='fake pedestal file')
parser.add_argument('-c', '--camera', type=int, help='camera.n')
parser.add_argument('-p', '--path', type=str, help='the base path to the output directory')
parser.add_argument('-n', '--nModules', type=int, help='nModules')
parser.add_argument('-d', '--detType', type=str, default='', help='Epix100, Epix10ka, Epix10kaQuad, Epix10ka2M, ...')
parser.add_argument('--maxNevents', type=int, default='666666', help='max number of events to analyze')
parser.add_argument('--skipNevents', type=int, default=0, help='max number of events to skip at the start of each step')
parser.add_argument('--configScript', type=str, default='experimentSuiteConfig.py', help='name of python config file to load if any')
parser.add_argument('--detObj', help='"raw", "calib", "image"')
parser.add_argument('-f','--file', type=str, help='run analysis only on file')
parser.add_argument('-L','--label', type=str, help='analysis label')
parser.add_argument('-t', '--threshold', help="threshold (ADU or keV or wave8) depending on --detObj")
parser.add_argument('--fluxCut', type=float, help="minimum flux to be included in analysis")
parser.add_argument('--special', type=str, help='comma-separated list of special behaviors - maybe this is too lazy. E.g. positiveParity,doKazEvents,...')
args = parser.parse_args()

##mymodule = importlib.import_module(full_module_name)

## for standalone analysis
Expand Down Expand Up @@ -153,27 +147,27 @@ def __init__(self, analysisType="scan"):
if args.fluxCut is not None:
self.fluxCut = args.fluxCut
try:
self.runRange = eval(args.runRange) ## in case needed
except Exception:
self.runRange = eval(args.runRange) ## in case needed
except:
self.runRange = None

self.fivePedestalRun = args.fivePedestalRun ## in case needed
self.fakePedestal = args.fakePedestal ## in case needed
self.fivePedestalRun = args.fivePedestalRun ## in case needed
self.fakePedestal = args.fakePedestal ## in case needed
if self.fakePedestal is not None:
self.fakePedestalFrame = np.load(self.fakePedestal) ##cast to uint32???

if args.detType == "":
self.fakePedestalFrame = np.load(self.fakePedestal) ##cast to uint32???
if args.detType == '':
## assume epix10k for now
if args.nModules is not None:
self.detType = self.ePix10k_cameraTypes[args.nModules]
else:
self.detType = args.detType

self.special = args.special
## done with configuration
## done with configuration

self.ds = None
self.det = None ## do we need multiple dets in an array? or self.secondDet?
self.det = None ## do we need multiple dets in an array? or self.secondDet?

##self.setupPsana()
##do this later or skip for -file
Expand All @@ -191,6 +185,7 @@ def setROI(self, roiFile=None, roi=None):
np.save(roiFile, roi)
self.ROI = roi


def noCommonModeCorrection(self, frame):
return frame

Expand Down Expand Up @@ -249,40 +244,31 @@ def colCommonModeCorrection(self, frame, arbitraryCut=1000):
rowOffset += self.detRowsPerBank
return frame

def commonModeCorrection(self, frame, arbitraryCut=1000):
## this takes a 2d frame
## cut keeps photons in common mode - e.g. set to <<1 photon

##rand = np.random.random()
for r in range(self.detRows):
colOffset = 0
##for b in range(0, self.detNbanks):
for b in range(0, 2):
try:
rowCM = np.median(
frame[r, colOffset : colOffset + self.detColsPerBank][
frame[r, colOffset : colOffset + self.detColsPerBank] < arbitraryCut
]
)
##if r == 280 and rand > 0.999:
##print(b, frame[r, colOffset:colOffset + self.detColsPerBank], \
# rowCM, rowCM<arbitraryCut-1, rowCM*(rowCM<arbitraryCut-1))
##frame[r, colOffset:colOffset + self.detColsPerBank] -= rowCM*(rowCM<arbitraryCut-1)
frame[r, colOffset : colOffset + self.detColsPerBank] -= rowCM
##if r == 280 and rand > 0.999:
##print(frame[r, colOffset:colOffset + self.detColsPerBank], \
# np.median(frame[r, colOffset:colOffset + self.detColsPerBank]))
except Exception:
rowCM = -666
print("rowCM problem")
print(frame[r, colOffset : colOffset + self.detColsPerBank])
colOffset += self.detColsPerBank
return frame
def isBeamEvent(self, evt):
ec = self.getEventCodes(evt)
##print(ec[280], ec[281], ec[282], ec[283], ec[284], ec[285] )
if ec[self.runCode]:
self.nRunCodeEvents += 1
if ec[self.daqCode]:
self.nDaqCodeEvents += 1
if self.fakeBeamCode:
return True
if ec[self.beamCode]:
self.nBeamCodeEvents += 1
return True
return False

def dumpEventCodeStatistics(self):
print("have counted %d run triggers, %d DAQ triggers, %d beam events"
%(self.nRunCodeEvents,
self.nDaqCodeEvents,
self.nBeamCodeEvents)
)

if __name__ == "__main__":
bSS = BasicSuiteScript()
print("have built a BasicSuiteScript")
bSS.setupPsana()
evt = bSS.getEvt()
print(dir(evt))

Loading

0 comments on commit 118bfea

Please sign in to comment.