Fixing the 4-component version of GST (#58)

* Correct the 4-component parameterization of GST

The original GST paper (http://doi.org/10.1007/s11467-013-0319-7) gives a parameterization of the CR flux in several populations based on p, He, C, O and Fe mass groups (Table 3). I see that this doesn't match the 5-component model that is used in some working groups (that's why I didn't touch the 5-component model) but for the 4-component version it would be appropriate to use the correct Z values and combine C and O.

* Update GST_IT

Updating the 4-component GST model to the 4-population version, having C and O combined with the correct Z and Fe combined with the super heavy component

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* Rename new flux model `GlobalFitGST4Comp` to keep _IT model the same

* Apply suggestion from @kjmeagher

---------

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
Co-authored-by: Kevin Meagher <11620178+kjmeagher@users.noreply.github.com>

Author: 3 people

src/simweights/__init__.py

@@ -28,6 +28,7 @@
     "GenieSurface",
     "GenieWeighter",
     "GlobalFitGST",
+    "GlobalFitGST4Comp",
     "GlobalFitGST_IT",
     "GlobalSplineFit",
     "GlobalSplineFit5Comp",
@@ -57,6 +58,7 @@
     GaisserH4a_IT,
     GaisserHillas,
     GlobalFitGST,
+    GlobalFitGST4Comp,
     GlobalFitGST_IT,
     GlobalSplineFit,
     GlobalSplineFit5Comp,

src/simweights/_fluxes.py

@@ -346,6 +346,39 @@ class GlobalFitGST_IT(CosmicRayFlux):
     )
 
 
+class GlobalFitGST4Comp(CosmicRayFlux):
+    r"""GlobalFitGST4Comp (4 populations) for four components [p, He, O, Fe].
+
+    The Oxygen group is the sum of Carbon and Oxygen groups.
+    The Iron group is the sum of Iron, Tellurium and Mercury groups.
+    Provided by J. Saffer.
+    """
+
+    pdgids = PDGID_4COMP
+    _funcs = (
+        lambda E: (
+            0.7000 * E**-2.66 * exp(-E / 1.2e5)
+            + 0.0150 * E**-2.4 * exp(-E / 4e6)
+            + 0.0012 * E**-2.4 * exp(-E / 1.5e9)
+            + 0.00012 * E**-2.4 * exp(-E / 40e9)
+        ),
+        lambda E: 0.3200 * E**-2.58 * exp(-E / 1.2e5 / 2) + 0.0065 * E**-2.3 * exp(-E / 4e6 / 2),
+        lambda E: (
+            0.0100 * E**-2.40 * exp(-E / 1.2e5 / 6)
+            + 0.0006 * E**-2.3 * exp(-E / 4e6 / 6)
+            + 0.0130 * E**-2.40 * exp(-E / 1.2e5 / 8)
+            + 0.0007 * E**-2.3 * exp(-E / 4e6 / 8)
+        ),
+        lambda E: (
+            0.0060 * E**-2.30 * exp(-E / 1.2e5 / 26)
+            + 0.00021 * E**-2.2 * exp(-E / 4e6 / 26)
+            + 0.0000011 * E**-2.2 * exp(-E / 1.5e9 / 26)
+            + 0.00001 * E**-2.2 * exp(-E / 4e6 / 52)
+            + 0.000053 * E**-2.2 * exp(-E / 4e6 / 80)
+        ),
+    )
+
+
 class GlobalSplineFitBase(CosmicRayFlux):
     r"""Data-driven spline fit of the cosmic ray spectrum by Dembinski et. al. \ [#GSFDembinski]_.