tpb: add alternative emission spectrum in PS matrix (#39)

src/legendoptics/data/tpb_polystrene_wlscomponent.dat

@@ -0,0 +1,98 @@
+# nm     dimensionless
+382.6768 0.0044
+385.7855 0.0263
+388.8944 0.0754
+392.0041 0.2049
+394.8323 0.4321
+397.0960 0.7382
+398.7947 1.0527
+400.2110 1.3921
+401.4858 1.7130
+402.6194 2.0374
+403.7534 2.3974
+404.8874 2.7611
+406.1627 3.1271
+407.2967 3.4955
+408.4313 3.9257
+409.5657 4.3250
+410.6592 4.6842
+411.5508 5.0309
+412.4017 5.3446
+413.2526 5.6583
+413.9905 5.9668
+415.0964 6.3547
+416.2303 6.7041
+417.3642 7.0642
+418.4980 7.4100
+419.6319 7.7594
+420.7658 8.1159
+421.8995 8.4581
+423.0330 8.7718
+424.1664 9.0784
+425.2998 9.3814
+426.5746 9.6944
+427.9906 10.0053
+429.5482 10.3380
+431.3883 10.6714
+433.5110 11.0034
+436.1979 11.2340
+439.3077 11.3687
+442.4167 11.4282
+445.5255 11.4721
+448.6336 11.4356
+451.7414 11.3680
+454.8489 11.2771
+457.9556 11.1084
+461.0618 10.8891
+464.0262 10.6095
+466.7076 10.2948
+469.1063 9.9674
+471.2226 9.6606
+473.1977 9.3651
+475.0313 9.0519
+476.7235 8.7190
+478.4158 8.3980
+480.2494 8.0790
+482.0830 7.7560
+483.9165 7.4350
+485.8912 7.0926
+487.8660 6.7603
+489.8407 6.4302
+491.8157 6.1142
+493.7906 5.8003
+495.9069 5.4893
+498.1643 5.1682
+500.2806 4.8587
+502.3969 4.5494
+504.6545 4.2444
+507.1946 3.9265
+509.8759 3.6080
+512.5573 3.2930
+515.5214 2.9898
+518.6272 2.7290
+521.7329 2.4527
+524.8387 2.1945
+527.9450 1.9817
+531.0515 1.7844
+534.1578 1.5781
+537.2644 1.3886
+540.3712 1.2354
+543.4783 1.0991
+546.5855 0.9732
+549.6930 0.8784
+552.8003 0.7706
+555.9078 0.6771
+559.0155 0.6030
+562.1233 0.5419
+565.2311 0.4821
+568.3390 0.4274
+571.4469 0.3780
+574.5549 0.3325
+577.6629 0.2934
+580.7710 0.2608
+583.8791 0.2269
+586.9872 0.1956
+590.0954 0.1721
+593.2037 0.1551
+596.3119 0.1381
+598.9963 0.1138

src/legendoptics/tpb.py

@@ -69,6 +69,19 @@ def tpb_wls_emission() -> tuple[Quantity, Quantity]:
     return readdatafile("tpb_vm2000_wlscomponent.dat")
 
 
+def tpb_polystrene_wls_emission() -> tuple[Quantity, Quantity]:
+    """WLS Emission spectrum for TPB in .
+
+    [Francini2013]_ measure the emission spectrum of TPB in a polystrene matrix
+    at an excitation wavelength of 128nm and at 87K, so exactly in our experimental
+    conditions. The major differences brougth by the embedding in the PS matrix is the
+    shift of the main emission peak and a loss of vibronic structures.
+
+    .. optics-plot::
+    """
+    return readdatafile("tpb_polystrene_wlscomponent.dat")
+
+
 def tpb_wls_absorption() -> tuple[Quantity, Quantity]:
     """Values reported in [Benson2018]_ for TPB evaporated on utraviolet-transmitting acrylic substrate.
 
@@ -92,7 +105,12 @@ def pyg4_tpb_attach_rindex(mat, reg) -> None:
         mat.addVecPropertyPint("RINDEX", λ.to("eV"), r)
 
 
-def pyg4_tpb_attach_wls(mat, reg, quantum_efficiency: bool | float = True) -> None:
+def pyg4_tpb_attach_wls(
+    mat,
+    reg,
+    quantum_efficiency: bool | float = True,
+    emission_spectrum: str = "default",
+) -> None:
     """Attach wavelength shifting properties to the given tpb material instance.
 
     Parameters
@@ -101,20 +119,30 @@ def pyg4_tpb_attach_wls(mat, reg, quantum_efficiency: bool | float = True) -> No
         If `False`, disable attaching any photon number information. If `True`, use the values
         from .tpb_quantum_efficiency. If specified as a number, directly attach this number as
         mean number of emitted photons.
+    emission_spectrum
+        either `default` or `polystrene_matrix`
 
     See Also
     --------
     .tpb_wls_absorption
     .tpb_wls_emission
+    .tpb_polystrene_wls_emission
     .tpb_wls_timeconstant
     .tpb_quantum_efficiency
     """
     from legendoptics.pyg4utils import pyg4_sample_λ
 
+    if emission_spectrum not in ["default", "polystrene_matrix"]:
+        raise ValueError("invalid parameter value of emission_spectrum")
+
+    emission_fn = tpb_wls_emission
+    if emission_spectrum == "polystrene_matrix":
+        emission_fn = tpb_polystrene_wls_emission
+
     λ_full = pyg4_sample_λ(112 * u.nm, 650 * u.nm, 800)
 
     absorption = InterpolatingGraph(*tpb_wls_absorption())(λ_full)
-    emission = InterpolatingGraph(*tpb_wls_emission())(λ_full)
+    emission = InterpolatingGraph(*emission_fn())(λ_full)
     # make sure that the scintillation spectrum is zero at the boundaries.
     emission[0] = 0
     emission[-1] = 0