Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add SAA, eclipse and flare flag in RHESSISummaryTimeSeries #124

Open
somusset opened this issue Jun 12, 2024 · 4 comments
Open

Add SAA, eclipse and flare flag in RHESSISummaryTimeSeries #124

somusset opened this issue Jun 12, 2024 · 4 comments

Comments

@somusset
Copy link

Describe the feature

It would be useful, when accessing the count rates from the RHESSI Summary observation files, to also get the most useful flags (to my opinion those are the SAA, eclipse and flare flags).

Proposed solution

As far as I understood the reading of this info will be done in the function parse_observing_summary_hdulist and this is how I accessed it (by imitation of the way the countrates are read):

    # calculate time array
    flag_reference_time_ut = parse_time(hdulist[14].data.field("UT_REF")[0], format="utime")
    flag_time_interval_sec = hdulist[14].data.field("TIME_INTV")[0]
    dim = hdulist[14].data.field("N_TIME_INTV")[0]

    flag_time_array = parse_time(flag_reference_time_ut) + TimeDelta(
        flag_time_interval_sec * np.arange(dim) * u.second)

    # extract flag arrays for SAA, ECLIPSE, FLARE
    flag_names = hdulist[14].data.field("FLAG_IDS")[0] 
    
    saa_flag = hdulist[15].data.field("FLAGS")[:,np.where(flag_names == "SAA_FLAG")[0][0]]
    eclipse_flag = hdulist[15].data.field("FLAGS")[:,np.where(flag_names == "ECLIPSE_FLAG")[0][0]]
    flare_flag = hdulist[15].data.field("FLAGS")[:,np.where(flag_names == "FLARE_FLAG")[0][0]]

One thing to note is that the time is redefined for the flags - I imagine that it should be the same time array as the one defined for the countrates, but since it is defined again for the flag, I guess that should be checked. Or maybe providing the two time arrays is not a problem.

@Abinash-bit
Copy link
Contributor

Abinash-bit commented Nov 29, 2024

def parse_observing_summary_hdulist(hdulist):
    """
    Parse a RHESSI observation summary file.

    Parameters
    ----------
    hdulist : `list`
        The HDU list from the fits file.

    Returns
    -------
    out : `dict`
        Returns a dictionary containing header and data information.
    """
    header = hdulist[0].header

    # Parse count rates
    reference_time_ut = parse_time(hdulist[5].data.field("UT_REF")[0], format="utime")
    time_interval_sec = hdulist[5].data.field("TIME_INTV")[0]
    dim = hdulist[5].data.field("N_TIME_INTV")[0]

    time_array = parse_time(reference_time_ut) + TimeDelta(
        time_interval_sec * np.arange(dim) * u.second
    )

    # Extract count rates
    compressed_countrate = np.array(hdulist[6].data.field("countrate"))
    countrate = uncompress_countrate(compressed_countrate)

    # Parse flag data
    flag_reference_time_ut = parse_time(hdulist[14].data.field("UT_REF")[0], format="utime")
    flag_time_interval_sec = hdulist[14].data.field("TIME_INTV")[0]
    flag_dim = hdulist[14].data.field("N_TIME_INTV")[0]

    flag_time_array = parse_time(flag_reference_time_ut) + TimeDelta(
        flag_time_interval_sec * np.arange(flag_dim) * u.second
    )

    # Extract flag arrays
    flag_names = hdulist[14].data.field("FLAG_IDS")[0]
    saa_flag = hdulist[15].data.field("FLAGS")[:, np.where(flag_names == "SAA_FLAG")[0][0]]
    eclipse_flag = hdulist[15].data.field("FLAGS")[:, np.where(flag_names == "ECLIPSE_FLAG")[0][0]]
    flare_flag = hdulist[15].data.field("FLAGS")[:, np.where(flag_names == "FLARE_FLAG")[0][0]]

    # Check if time arrays match
    if not np.array_equal(time_array, flag_time_array):
        Warning.warn("Time arrays for count rates and flags do not match. Proceed with caution.")

    # Organize data
    data = {
        "time": time_array,
        "countrate": countrate,
        "labels": [
            "3 - 6 keV",
            "6 - 12 keV",
            "12 - 25 keV",
            "25 - 50 keV",
            "50 - 100 keV",
            "100 - 300 keV",
            "300 - 800 keV",
            "800 - 7000 keV",
            "7000 - 20000 keV",
        ],
        "flags": {
            "time": flag_time_array,
            "saa_flag": saa_flag,
            "eclipse_flag": eclipse_flag,
            "flare_flag": flare_flag,
        },
    }

    return header, data

I have changed the parse_observing_summary_hdulist in rhessi.py but i can not visualize the output can anyone help.

@nabobalis
Copy link
Contributor

Hi @Abinash-bit, what do you mean by visualize the output?

@ianan
Copy link

ianan commented Dec 18, 2024

Some thoughts

  • I'm sure the flag binning is the same as the count rates, so don't need the additional time bins
  • Would also be useful to include the attenuator flag (which would be have values 0, 1, 3) its ATTENUATOR_STATE or hdulist[15].data.field("flags")[14]
  • I'm not sure the flare flags in these files actually contain any flare time, I've checked a few files when I know there are flares and the flags are all 0.... am investigating....update: the flare flag in these files are all 0, so would need to get that info from somewhere else.
  • In the above code I'd maybe change labels to just use the info in the files, i.e. labels=[ff+' keV' for ff in hdulist[5].data[0][7]]

@ianan
Copy link

ianan commented Dec 20, 2024

Further, that the flare flags appearing in the RHESSI quicklook plots are taken from the RHESSI flare list and not from the flare info in these summary files, as what is called the flare flare is just all 0. Specifically in sswidl hsi_obs_summary() gets the flags via his_obs_summ_flag__define https://hesperia.gsfc.nasa.gov/ssw/hessi/idl/qlook_archive/hsi_flare_list__define.pro getting some flags from the summary files but uses the reset_flare_flag method which is reading the flare catalogue to determine when the flare flag is changing.

So to have this functionality in python need to read the flare list, search and return particular flag info - it seems this hasn't been implemented yet as per #36 ?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

4 participants