Use WCSim properties for water with the lbne detector

1 files changed, 74 insertions, 0 deletions

diff --git a/optics.py b/optics.py
index d65c5f9..c31d209 100644
--- a/optics.py
+++ b/optics.py
@@ -47,6 +47,8 @@ r7081hqe_photocathode.absorb = r7081hqe_photocathode.detect
 # remaining photons are diffusely reflected
 r7081hqe_photocathode.set('reflect_diffuse', 1.0 - r7081hqe_photocathode.detect[:,1] - r7081hqe_photocathode.absorb[:,1], wavelengths=r7081hqe_photocathode.detect[:,0])
 
+######################## SNO+ materials ##############################
+
 # water data comes from 'lightwater_sno' material in the SNO+ optics database
 water = Material('water')
 water.density = 1.0 # g/cm^3
@@ -168,3 +170,75 @@ labppo_scintillator.set('absorption_length', wavelengths=[ 200.0, 215.0, 225.0,
 labppo_scintillator.set('absorption_wls_length', wavelengths=[ 200.0, 215.0, 225.0, 235.0, 245.0, 255.0, 265.0, 275.0, 285.0, 295.0, 305.0, 315.0, 325.0, 335.0, 345.0, 355.0, 365.0, 375.0, 385.0, 395.0, 405.0, 415.0, 425.0, 435.0, 445.0, 455.0, 465.0, 475.0, 485.0, 495.0, 505.0, 515.0, 525.0, 535.0, 545.0, 555.0, 565.0, 575.0, 585.0, 595.0, 800.0 ], value=[ 0.00109808, 0.000693615, 0.000437277, 0.000276096, 0.000173933, 0.000103044, 6.03065e-05, 3.56199e-05, 2.25651e-05, 1.60844e-05, 1.62228e-05, 2.21832e-05, 3.03067e-05, 5.24247e-05, 0.000484492, 0.0127541, 0.294651, 4.15892, 15.6626, 25.4821, 36.1763, 50.723, 63.9962, 77.7386, 89.8611, 103.551, 115.859, 125.442, 146.664, 163.748, 194.986, 194.652, 202.262, 234.033, 262.861, 279.883, 323.991, 299.632, 287.367, 303.42, 303.42 ])
 labppo_scintillator.set('scattering_length', wavelengths=[ 200.0, 215.0, 225.0, 235.0, 245.0, 255.0, 265.0, 275.0, 285.0, 295.0, 305.0, 315.0, 325.0, 335.0, 345.0, 355.0, 365.0, 375.0, 385.0, 395.0, 405.0, 415.0, 425.0, 435.0, 445.0, 455.0, 465.0, 475.0, 485.0, 495.0, 505.0, 515.0, 525.0, 535.0, 545.0, 555.0, 565.0, 575.0, 585.0, 595.0, 800.0 ], value=[ 0.290404, 0.387993, 0.465472, 0.553999, 0.654579, 0.768259, 0.89613, 1.03933, 1.19903, 1.37645, 1.57287, 1.78958, 2.02794, 2.28936, 2.57527, 2.88715, 3.22652, 3.59497, 3.9941, 4.42558, 4.89108, 5.3924, 5.93129, 6.50958, 7.12918, 7.79198, 8.49995, 9.2551, 10.0595, 10.9152, 11.8243, 12.7892, 13.8118, 14.8947, 16.04, 17.25, 18.5273, 19.8742, 21.2933, 22.7871, 22.7871 ])
 
+################### WCSim materials #####################
+
+water_wcsim = Material('water')
+water_wcsim.density = 1.0 # g/cm^3
+water_wcsim.composition = { 'H' : 0.1119, 'O' : 0.8881 } # fraction by mass
+hc_over_GeV = 1.2398424468024265e-06 # h_Planck * c_light / GeV / nanometer
+wcsim_wavelengths = hc_over_GeV / np.array([ 1.56962e-09, 1.58974e-09, 1.61039e-09, 1.63157e-09, 
+       1.65333e-09, 1.67567e-09, 1.69863e-09, 1.72222e-09, 
+       1.74647e-09, 1.77142e-09,1.7971e-09, 1.82352e-09, 
+       1.85074e-09, 1.87878e-09, 1.90769e-09, 1.93749e-09, 
+       1.96825e-09, 1.99999e-09, 2.03278e-09, 2.06666e-09,
+       2.10169e-09, 2.13793e-09, 2.17543e-09, 2.21428e-09, 
+       2.25454e-09, 2.29629e-09, 2.33962e-09, 2.38461e-09, 
+       2.43137e-09, 2.47999e-09, 2.53061e-09, 2.58333e-09, 
+       2.63829e-09, 2.69565e-09, 2.75555e-09, 2.81817e-09, 
+       2.88371e-09, 2.95237e-09, 3.02438e-09, 3.09999e-09,
+       3.17948e-09, 3.26315e-09, 3.35134e-09, 3.44444e-09, 
+       3.54285e-09, 3.64705e-09, 3.75757e-09, 3.87499e-09, 
+       3.99999e-09, 4.13332e-09, 4.27585e-09, 4.42856e-09, 
+       4.59258e-09, 4.76922e-09, 4.95999e-09, 5.16665e-09, 
+       5.39129e-09, 5.63635e-09, 5.90475e-09, 6.19998e-09 ])[::-1] #reversed
+
+water_wcsim.set('refractive_index', 
+                wavelengths=wcsim_wavelengths,
+                value=np.array([1.32885, 1.32906, 1.32927, 1.32948, 1.3297, 1.32992, 1.33014, 
+                          1.33037, 1.3306, 1.33084, 1.33109, 1.33134, 1.3316, 1.33186, 1.33213,
+                          1.33241, 1.3327, 1.33299, 1.33329, 1.33361, 1.33393, 1.33427, 1.33462,
+                          1.33498, 1.33536, 1.33576, 1.33617, 1.3366, 1.33705, 1.33753, 1.33803,
+                          1.33855, 1.33911, 1.3397, 1.34033, 1.341, 1.34172, 1.34248, 1.34331,
+                          1.34419, 1.34515, 1.3462, 1.34733, 1.34858, 1.34994, 1.35145, 1.35312,
+                          1.35498, 1.35707, 1.35943, 1.36211, 1.36518, 1.36872, 1.37287, 1.37776,
+                          1.38362, 1.39074, 1.39956, 1.41075, 1.42535])[::-1] #reversed
+)
+water_wcsim.set('absorption_length',
+                wavelengths=wcsim_wavelengths,
+                value=np.array([22.8154, 28.6144, 35.9923, 45.4086, 57.4650,
+                 72.9526, 75,      81.2317, 120.901, 160.243,
+                 193.797, 215.045, 227.786, 243.893, 294.113,
+                 321.735, 342.931, 362.967, 378.212, 449.602,
+                 740.143, 1116.06, 1438.78, 1615.48, 1769.86,
+                 2109.67, 2304.13, 2444.97, 3076.83, 4901.5,
+                 6666.57, 7873.95, 9433.81, 10214.5, 10845.8,
+                 15746.9, 20201.8, 22025.8, 21142.2, 15083.9,
+                 11751,   8795.34, 8741.23, 7102.37, 6060.68,
+                 4498.56, 3039.56, 2232.2,  1938,    1811.58,
+                 1610.32, 1338.7,  1095.3,  977.525, 965.258,
+                 1082.86, 876.434, 633.723, 389.87,  142.011])[::-1] / 100.0 # reversed, cm->m
+                )
+      
+water_wcsim.set('scattering_length',
+                wavelengths=wcsim_wavelengths,
+                value=np.array([167024.4, 158726.7, 150742,
+                          143062.5, 135680.2, 128587.4,
+                          121776.3, 115239.5, 108969.5,
+                          102958.8, 97200.35, 91686.86,
+                          86411.33, 81366.79, 76546.42,
+                          71943.46, 67551.29, 63363.36,
+                          59373.25, 55574.61, 51961.24,
+                          48527.00, 45265.87, 42171.94,
+                          39239.39, 36462.50, 33835.68,
+                          31353.41, 29010.30, 26801.03,
+                          24720.42, 22763.36, 20924.88,
+                          19200.07, 17584.16, 16072.45,
+                          14660.38, 13343.46, 12117.33,
+                          10977.70, 9920.416, 8941.407,
+                          8036.711, 7202.470, 6434.927,
+                          5730.429, 5085.425, 4496.467,
+                          3960.210, 3473.413, 3032.937,
+                          2635.746, 2278.907, 1959.588,
+                          1675.064, 1422.710, 1200.004,
+                          1004.528, 833.9666, 686.1063])[::-1] / 100.0 * 0.625 # reversed, cm -> m, * magic tuning constant
+                )