Skip to content

Commit

Permalink
oh, weird inheritance didn't work?
Browse files Browse the repository at this point in the history
  • Loading branch information
josh0-jrg committed Sep 20, 2024
1 parent 8b3169f commit 3fb0932
Show file tree
Hide file tree
Showing 2 changed files with 77 additions and 77 deletions.
39 changes: 3 additions & 36 deletions flamedisx/lz/lz.py
Original file line number Diff line number Diff line change
Expand Up @@ -64,8 +64,8 @@ def build_position_map_from_data(map_file, axis_names, bins):


class LZSource:
path_s1_corr_LZAP = 'new_data/s1Area_Correction_TPC_SR3_06Apr23.json'
path_s2_corr_LZAP = 'new_data/s2Area_Correction_TPC_SR3_06Apr23.json'
path_s1_corr_LZAP = 'new_data/s1Area_Correction_TPC_SR3_radon_31Jan2024.json'#'new_data/s1Area_Correction_TPC_SR3_06Apr23.json'
path_s2_corr_LZAP = 'new_data/s2Area_Correction_TPC_SR3_radon_31Jan2024.json'#'new_data/s2Area_Correction_TPC_SR3_06Apr23.json'
path_s1_corr_latest = 'new_data/s1Area_Correction_TPC_SR3_radon_31Jan2024.json'
path_s2_corr_latest = 'new_data/s2Area_Correction_TPC_SR3_radon_31Jan2024.json'

Expand Down Expand Up @@ -340,40 +340,7 @@ def __init__(self, *args, **kwargs):
if ('detector' not in kwargs):
kwargs['detector'] = 'lz_SR3'
super().__init__(*args, **kwargs)
def mean_yield_electron(self, energy):
#Refactored to take constants from lzlama !397
m1=12.4886
m2=85.0
m3=0.6050
m4= 2.14687
m5=25.721
m6=-1.0
m7=59.651,
m8=3.6869
m9=0.2872
m10=0.1121
Wq_eV = self.Wq_keV * 1e3

Nq = energy * 1e3 / Wq_eV


Qy = m1 + (m2 - m1) / pow((1. + pow(energy /m3,m4)),m9) + \
m5 + (m6 - m5) / pow((1. + pow(energy /m7, m8)), m10)

coeff_TI = tf.cast(pow(1. / XENON_REF_DENSITY, 0.3), fd.float_type())
coeff_Ni = tf.cast(pow(1. / XENON_REF_DENSITY, 1.4), fd.float_type())
coeff_OL = tf.cast(pow(1. / XENON_REF_DENSITY, -1.7) /
fd.tf_log10(1. + coeff_TI * coeff_Ni * pow(XENON_REF_DENSITY, 1.7)), fd.float_type())

Qy *= coeff_OL * fd.tf_log10(1. + coeff_TI * coeff_Ni * pow(self.density, 1.7)) * pow(self.density, -1.7)

nel_temp = Qy * energy
# Don't let number of electrons go negative
nel = tf.where(nel_temp < 0,
0 * nel_temp,
nel_temp)

return nel


@export
class LZGammaSource(LZSource, fd.nest.nestGammaSource):
Expand Down
115 changes: 74 additions & 41 deletions flamedisx/nest/lxe_sources.py
Original file line number Diff line number Diff line change
Expand Up @@ -217,52 +217,25 @@ def __init__(self, *args, energy_min=0.01, energy_max=10., num_energies=1000, en
fd_nest.MakeS2)

# quanta_splitting.py

def mean_yield_electron(self, energy):
#NEED TO REFACTOR IN TERMS OF M1,M2...,M10 ignoring functional forms
er_m1_a=30.66
er_m1_b=6.1978
er_m1_d=73.855
er_m1_e=2.0318
er_m1_f=0.41883
er_m10_a=0.0508273937
er_m10_b=0.1166087199
er_m10_c= 0.0508273937
er_m10_d=1.39260460e+02
er_m10_e=-0.65763592
er_Qy_a=77.2931084
er_Qy_b=0.13946236
er_Qy_c=0.52561312
er_Qy_d=1.82217496
er_Qy_e=2.82528809
er_Qy_f=1.82217496
er_Qy_g=144.65029656
er_Qy_h=-2.80532006
er_Qy_i=0.3344049589
er_Qy_k=7.02921301
er_Qy_l=98.27936794
er_Qy_m=7.0292130
er_Qy_n=256.48156448
er_Qy_o=1.29119251
er_Qy_p=4.285781736
#Refactored to take constants from lzlama !397
m1=12.4886
m2=85.0
m3=0.6050
m4= 2.14687
m5=25.721
m6=-1.0
m7=59.651
m8=3.6869
m9=0.2872
m10=0.1121
Wq_eV = self.Wq_keV * 1e3

Nq = energy * 1e3 / Wq_eV

m1 = tf.cast(er_m1_a + (er_m1_b - er_m1_a) / pow(1. + pow(self.drift_field / er_m1_d, er_m1_e), er_m1_f),
fd.float_type())
m5 = tf.cast(Nq / energy / (1 + self.alpha * tf.math.erf(0.05 * energy)), fd.float_type()) - m1
m10 = tf.cast((er_m10_a + (er_m10_b - er_m10_c) /
(1 + pow(self.drift_field / er_m10_d, er_m10_e))),
fd.float_type())

Qy = m1 + (er_Qy_a - m1) / pow((1. + pow(energy / (fd.tf_log10(tf.cast(self.drift_field, fd.float_type())) *
er_Qy_b + er_Qy_c),
er_Qy_d + (er_Qy_e - er_Qy_f) /
(1 + pow(self.drift_field / er_Qy_g, er_Qy_h)))),
er_Qy_i) + \
m5 + (0. - m5) / pow((1. + pow(energy / (er_Qy_k + (er_Qy_l - er_Qy_m) /
(1. + pow(self.drift_field / er_Qy_n, er_Qy_o))), er_Qy_p)), m10)

Qy = m1 + (m2 - m1) / pow((1. + pow(energy /m3,m4)),m9) + \
m5 + (m6 - m5) / pow((1. + pow(energy /m7, m8)), m10)

coeff_TI = tf.cast(pow(1. / XENON_REF_DENSITY, 0.3), fd.float_type())
coeff_Ni = tf.cast(pow(1. / XENON_REF_DENSITY, 1.4), fd.float_type())
Expand All @@ -278,6 +251,66 @@ def mean_yield_electron(self, energy):
nel_temp)

return nel
# def mean_yield_electron(self, energy):
# #NEED TO REFACTOR IN TERMS OF M1,M2...,M10 ignoring functional forms
# er_m1_a=30.66
# er_m1_b=6.1978
# er_m1_d=73.855
# er_m1_e=2.0318
# er_m1_f=0.41883
# er_m10_a=0.0508273937
# er_m10_b=0.1166087199
# er_m10_c= 0.0508273937
# er_m10_d=1.39260460e+02
# er_m10_e=-0.65763592
# er_Qy_a=77.2931084
# er_Qy_b=0.13946236
# er_Qy_c=0.52561312
# er_Qy_d=1.82217496
# er_Qy_e=2.82528809
# er_Qy_f=1.82217496
# er_Qy_g=144.65029656
# er_Qy_h=-2.80532006
# er_Qy_i=0.3344049589
# er_Qy_k=7.02921301
# er_Qy_l=98.27936794
# er_Qy_m=7.0292130
# er_Qy_n=256.48156448
# er_Qy_o=1.29119251
# er_Qy_p=4.285781736
# Wq_eV = self.Wq_keV * 1e3

# Nq = energy * 1e3 / Wq_eV

# m1 = tf.cast(er_m1_a + (er_m1_b - er_m1_a) / pow(1. + pow(self.drift_field / er_m1_d, er_m1_e), er_m1_f),
# fd.float_type())
# m5 = tf.cast(Nq / energy / (1 + self.alpha * tf.math.erf(0.05 * energy)), fd.float_type()) - m1
# m10 = tf.cast((er_m10_a + (er_m10_b - er_m10_c) /
# (1 + pow(self.drift_field / er_m10_d, er_m10_e))),
# fd.float_type())

# Qy = m1 + (er_Qy_a - m1) / pow((1. + pow(energy / (fd.tf_log10(tf.cast(self.drift_field, fd.float_type())) *
# er_Qy_b + er_Qy_c),
# er_Qy_d + (er_Qy_e - er_Qy_f) /
# (1 + pow(self.drift_field / er_Qy_g, er_Qy_h)))),
# er_Qy_i) + \
# m5 + (0. - m5) / pow((1. + pow(energy / (er_Qy_k + (er_Qy_l - er_Qy_m) /
# (1. + pow(self.drift_field / er_Qy_n, er_Qy_o))), er_Qy_p)), m10)

# coeff_TI = tf.cast(pow(1. / XENON_REF_DENSITY, 0.3), fd.float_type())
# coeff_Ni = tf.cast(pow(1. / XENON_REF_DENSITY, 1.4), fd.float_type())
# coeff_OL = tf.cast(pow(1. / XENON_REF_DENSITY, -1.7) /
# fd.tf_log10(1. + coeff_TI * coeff_Ni * pow(XENON_REF_DENSITY, 1.7)), fd.float_type())

# Qy *= coeff_OL * fd.tf_log10(1. + coeff_TI * coeff_Ni * pow(self.density, 1.7)) * pow(self.density, -1.7)

# nel_temp = Qy * energy
# # Don't let number of electrons go negative
# nel = tf.where(nel_temp < 0,
# 0 * nel_temp,
# nel_temp)

# return nel

def mean_yield_quanta(self, *args):
energy = args[0]
Expand Down

0 comments on commit 3fb0932

Please sign in to comment.