Handling the SAA passages

[GallegoSav]

For DC3 the SAA passages are simulated, so for sources and background we need to remove events when timestamp correspond to a SAA passage.

I was thinking that @ckarwin could do this cut during the conversion of the tra file into fits file. I share a example below using pandas which is quiet convenient for this type of tasks.

import pandas as pd

#we read the SAA LC
SAAProtonLC = pd.read_csv("/lustre/project/nhr-cosi/MC/DC3/INPUT/LightCurve/SAAproton_4MeVto2000.dat",skiprows=1,engine="python",skipfooter=2,names=["DP","time","flux"],sep=" ")
SAAProtonLC.set_index("time",inplace=True)


#your background or source file you read into a dataframe
df = pd.DataFrame(....)
df.set_index("timestamp",inplace=True)


#now we reindex the SAALC file to the dataframe time
SAAProtonLC_tmp = SAAProtonLC.reindex(df.index,method="nearest")
df["SAAflux"] = SAAProtonLC_tmp["flux"]

#we remove all the events where timestamp correspond to SAA flux >0
df = df[df["SAAflux"]==0]

Of course the exposure time will need to be corrected (for 3 months you remove ~16 days) but you can compute this once so we could create one function that compute the correct exposure time for a specific time interval. Here a example

df = pd.read_csv("/localscratch/sgallego/linkToXauron/MC/Background/COSIMEX/DC3/INPUT/LightCurves/SAAproton_4MeVto2000.dat",skiprows=1,engine="python",skipinitialspace=True,skipfooter=2,names=["DP","time","flux"],sep=" ")

# time period you want to study, I commented so this case will be the 3 months
#df = df[  (df["time"]<=TimeMin)  & (df["time"]<=TimeMax)    ]


df['datetime']= pd.to_datetime(df["time"],unit="s")
df.set_index("datetime",inplace=True)

# Calculate time deltas for each row
time_deltas = df.index.to_series().diff().fillna(pd.Timedelta(seconds=0))

#print(time_deltas)

# Identify periods where SAA_flux > 0 and where SAA_flux == 0
exposure_active = time_deltas[df['flux'] > 0].sum()  # Total exposure time for SAA_flux > 0
exposure_inactive = time_deltas[df['flux'] == 0].sum()  # Total exposure time for SAA_flux == 0

# Convert the Timedelta result to seconds, minutes, or hours as needed
exposure_active_seconds = exposure_active.total_seconds()
exposure_inactive_seconds = exposure_inactive.total_seconds()

print("Exposure time (SAA flux > 0):", exposure_active)
print("Exposure time (SAA flux == 0):", exposure_inactive)

print("Exposure time in seconds (SAA flux > 0):", exposure_active_seconds)
print("Exposure time in seconds (SAA flux == 0):", exposure_inactive_seconds)

Exposure time (SAA flux > 0): 16 days 05:00:00
Exposure time (SAA flux == 0): 76 days 03:39:00
Exposure time in seconds (SAA flux > 0): 1400400.0
Exposure time in seconds (SAA flux == 0): 6579540.0

[ckarwin]

Hi Savitri,

I think the cuts for the data can be handled with the dataIO class. We just need the time intervals for which COSI will be turned off during the SAA passage. This information should be part of the spacecraft file, and maybe for DC3 we can already add it to the ori file.

The change in exposure time can also be handled in cosipy. This can by done by setting the weights to zero during the scatt_map calculation, similar to what we do for the Earth occultation. Again, we just need the time intervals corresponding to the SAA passage.

I think @israelmcmc has some ideas for this. Let's all discuss this more at tomorrow's cosipy meeting.