Merge pull request #19 from jeremy-baier/sensitivity-evolution

Add sensitivity forecasting and detection probability calculations

.github/workflows/ci_tests.yml

@@ -17,13 +17,13 @@ jobs:
       fail-fast: false
       matrix:
         os: [ubuntu-latest, macos-latest]
-        python-version: [3.6, 3.7, 3.8, 3.9]
+        python-version: ["3.9", "3.10", "3.11", "3.12", "3.13"]
 
     steps:
     - name: Checkout repository
-      uses: actions/checkout@v2
+      uses: actions/checkout@v3
     - name: Set up Python ${{ matrix.python-version }}
-      uses: actions/setup-python@v2
+      uses: actions/setup-python@v4
       with:
         python-version: ${{ matrix.python-version }}
     - name: Install dependencies and package
@@ -53,11 +53,11 @@ jobs:
     # if: github.event_name == 'release'
 
     steps:
-    - uses: actions/checkout@v2
+    - uses: actions/checkout@v3
     - name: Set up Python
-      uses: actions/setup-python@v2
+      uses: actions/setup-python@v4
       with:
-        python-version: "3.9"
+        python-version: "3.10"
     - name: Build
       run: |
         python -m pip install --upgrade pip setuptools wheel
@@ -84,7 +84,7 @@ jobs:
         venv-wheel/bin/python -m pip install --upgrade pip setuptools
         venv-wheel/bin/python -m pip install ../dist/hasasia*.whl
         venv-wheel/bin/python -c "import hasasia;print(hasasia.__version__)"
-    - uses: actions/upload-artifact@v2
+    - uses: actions/upload-artifact@v4
       with:
         name: dist
         path: dist/*
@@ -94,17 +94,17 @@ jobs:
     runs-on: ubuntu-latest
     if: github.event_name == 'release'
     steps:
-    - uses: actions/checkout@v2
+    - uses: actions/checkout@v3
     - name: Set up Python
-      uses: actions/setup-python@v2
+      uses: actions/setup-python@v4
       with:
-        python-version: '3.9'
+        python-version: '3.10'
     - name: Install dependencies
       run: |
         python -m pip install --upgrade pip
         pip install setuptools wheel twine
     - name: Download wheel/dist from build
-      uses: actions/download-artifact@v2
+      uses: actions/download-artifact@v3
       with:
         name: dist
         path: dist

.gitignore

@@ -1,8 +1,12 @@
+#MacBook stuff
+.DS_Store
+.gitignore
+
 # Byte-compiled / optimized / DLL files
 __pycache__/
 *.py[cod]
 *$py.class
-
+**/.DS_Store
 # C extensions
 *.so
 

hasasia/forecast.py

@@ -0,0 +1,62 @@
+# -*- coding: utf-8 -*-
+import copy
+from hasasia import sensitivity as hsen
+
+yr_sec = 365.25*24*3600
+
+def get_sliced_spectra(psrs,
+                       freqs,
+                       start_mjd=0.0,
+                       end_mjd = 1_000_000,
+                       min_tspan_cut = 3,
+                       verbose=True,
+                       ):
+    """
+    Parameters 
+    ----------
+    
+    psrs : list of enterprise.Pulsar or list of hasasia.Pulsar objects
+        List of enterprise/hasasia Pulsar objects
+
+    freqs : array
+        Frequency array for the PTA in Hz
+
+    start_mjd : float
+        Start time in MJD for the slice
+
+    end_mjd : float
+        End time in MJD for the slice
+    
+    min_tspan_cut : float
+        Minimum time span in years for the pulsar to be included in the PTA spectrum
+    
+    verbose : bool
+        Print out the number of pulsars and the time span of the PTA spectrum
+
+    Return
+    ------
+    spectra : list of hasasia.Spectrum objects
+    """
+    psrs_post_cut = [] # list of pulsars after timespan cuts
+    psrs_copy = copy.deepcopy(psrs)
+    for _ , psr in enumerate(psrs_copy):
+        # filter the data around to the appropriate slice
+        psr.filter_data(start_time=start_mjd, end_time=end_mjd)
+        # If there are no TOAs remaining, or the time span of the pulsar is <3 cut.
+        if (psr.toas.size == 0):
+            pass # don't include pulsars without TOAs
+        elif min_tspan_cut is None:
+            psrs_post_cut.append(psr)
+        elif min_tspan_cut is not None and  (hsen.get_Tspan([psr]) < min_tspan_cut*yr_sec):
+            pass # don't include pulsars with less than min_tspan_cut years of data
+        else:
+            psrs_post_cut.append(psr)
+    spectra = []
+    for _ , p in enumerate(psrs_post_cut):
+        sp = hsen.Spectrum(p,freqs=freqs)
+        sp.name = p.name
+        _ = sp.NcalInv
+        spectra.append(sp)
+    if verbose is True:
+        print(f"PTA spectrum with {len(spectra)} psrs and Tspan {round(hsen.get_Tspan([psrs_post_cut])/yr_sec, 2)} created.")
+    return spectra