Skip to content

Adding features to the HDF file generation scripts #267

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.

Sign up for GitHub

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

Open
eirikgje wants to merge 5 commits into
base: akari
Choose a base branch
from
Open

Adding features to the HDF file generation scripts #267

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
f09bdae
Add solar-centric pixels to HDF files
eirikgje May 19, 2026
8ebd260
Enable tod and todz compression in the same file
eirikgje May 19, 2026
b5285b4
Enable extending the reset flag to 4 samples
eirikgje May 19, 2026
c92915e
Update fits2hdf to include previous changes
eirikgje May 28, 2026
cad5f41
Merge branch 'akari' into akari_hdf_updates
eirikgje May 28, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
66 changes: 36 additions & 30 deletions commander3/todscripts/AKARI/akari_commander_data_adapter.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -264,17 +264,13 @@ def fits2output_formatter(file, start_index, end_index, band,
# an okay assumption for now, but I just note it here.
start_det_ind = start_det_inds[band]
end_det_ind = start_det_ind + num_detectors[band]
if should_compress:
tot_adu = file[1]['DET'].read()[start_index:end_index, start_det_ind:end_det_ind]
else:
tot_flux = file[1]['FLUX'].read()[start_index:end_index, start_det_ind:end_det_ind]
tot_adu = file[1]['DET'].read()[start_index:end_index, start_det_ind:end_det_ind]
tot_flux = file[1]['FLUX'].read()[start_index:end_index, start_det_ind:end_det_ind]
tot_pixflag = file[1]['PIX_FLAG'].read()[start_index:end_index,
start_det_ind:end_det_ind, :].astype(bool)
for local_det_idx, detname in enumerate(band_dets[band]):
if should_compress:
out_data[f'{detname}/todz'] = tot_adu[:, local_det_idx]
else:
out_data[f'{detname}/tod'] = tot_flux[:, local_det_idx]
out_data[f'{detname}/tod'] = tot_flux[:, local_det_idx]
out_data[f'{detname}/todz'] = tot_adu[:, local_det_idx]
out_data[f'{detname}/pixel_flag'] = tot_pixflag[:, local_det_idx, :]
out_data['status_flag'] = file[1]['STATUS'].read()[start_index:end_index].astype(bool)
# NB! This is necessary because we interpret the shtop flag inversely internally in commander
Expand All @@ -300,11 +296,13 @@ class AKARICommanderDataAdapter(CommanderDataAdapter):
"""
def __init__(self, akari_fits_dir, nside, bands=BANDS,
num_detectors=NUM_DETECTORS, band_dets=BAND_DETS,
reference_time=REFERENCE_TIME):
reference_time=REFERENCE_TIME,
extend_reset_flag=False):
self.bands = bands
self.num_detectors = num_detectors
self.band_dets = band_dets
self.nside = nside
self.extend_reset_flag = extend_reset_flag


self.fsamp = {
Expand Down Expand Up @@ -341,8 +339,6 @@ def __init__(self, akari_fits_dir, nside, bands=BANDS,
self.nchunks[band] = len(self.chunk_file_map[band])
self.reference_time = reference_time

self.should_compress = True


def _calculate_chunk_files(self):
"""
Expand Down Expand Up @@ -410,8 +406,7 @@ def _process_akari_tod_reader_data(self):
# within a file rather than at the beginning of the file.
start_idx = self.todreaders[band].get_file_index_range(files[0])[0]
end_idx = self.todreaders[band].get_file_index_range(files[-1])[1]
curr_data = self.todreaders[band].get_data(start_idx, end_idx,
should_compress=self.should_compress)
curr_data = self.todreaders[band].get_data(start_idx, end_idx)
mode = curr_data['packet_id']
if len(files) > 1:
raise NotImplementedError("""We currently don't support loading gads flags for several
Expand All @@ -423,20 +418,23 @@ def _process_akari_tod_reader_data(self):
spin_axis_arr = np.zeros((2, len(self.band_dets[band])))
for detidx, det in enumerate(self.band_dets[band]):
todreader_data[band][det] = {}
if self.should_compress:
todreader_data[band][det]['todz'] = curr_data[f'{det}/todz']
else:
todreader_data[band][det]['tod'] = curr_data[f'{det}/tod']
todreader_data[band][det]['tod'] = curr_data[f'{det}/tod']
todreader_data[band][det]['todz'] = curr_data[f'{det}/todz']
todreader_data[band][det]['flag'] = self._process_akari_flags(
curr_data[f'{det}/pixel_flag'],
curr_data['frame_flag'],
curr_data['status_flag'], gads_flags[det], mode)
c = SkyCoord(ra=detlons[det]*u.deg, dec=detlats[det]*u.deg, frame='icrs')
curr_data['status_flag'], gads_flags[det], mode,
extend_reset_flag=self.extend_reset_flag)
c = SkyCoord(ra=detlons[det]*u.deg, dec=detlats[det]*u.deg, frame='icrs', distance=1*u.AU)
todreader_data[band][det]['pix'] = healpy.ang2pix(self.nside,
c.galactic.l.value,
c.galactic.b.value,
lonlat=True)

theta = 0.5 * np.pi - c.hcrs.dec.radian
phi = c.hcrs.ra.radian
todreader_data[band][det]['pix_solarcentric'] = healpy.ang2pix(
self.nside, theta, phi)
vecs = c.galactic.cartesian.get_xyz().value.T
spin_axis_arr[:,detidx] = self.todreaders[band]._ring_spin_axis(vecs)

Expand All @@ -448,15 +446,19 @@ def _process_akari_tod_reader_data(self):
mu = spin_axis_arr.mean(axis=1)
todreader_data[band]['spin_axis'] = mu


ra = curr_data['ra']
dec = curr_data['dec']
c = SkyCoord(ra=ra*u.deg, dec=dec*u.deg, frame='icrs')
c = SkyCoord(ra=ra*u.deg, dec=dec*u.deg, frame='icrs', distance=1*u.AU)
todreader_data[band]['pix'] = healpy.ang2pix(self.nside,
c.galactic.l.value,
c.galactic.b.value,
lonlat=True)

theta = 0.5 * np.pi - c.hcrs.dec.radian
phi = c.hcrs.ra.radian
todreader_data[band]['pix_solarcentric'] = healpy.ang2pix(
self.nside, theta, phi)



if 'start_time' not in todreader_data:
Expand All @@ -482,10 +484,11 @@ def _process_akari_tod_reader_data(self):
return todreader_data


def _process_akari_flags(self, pixflag_array, frameflag_array, statusflag_array, gadsflag_array,
def _process_akari_flags(self, pixflag_array, frameflag_array,
statusflag_array, gadsflag_array,
mode, pixflagmap=PIX_FLAG_MAP, frameflagmap=FRAME_FLAG_MAP,
statusflagmap=STATUS_FLAG_MAP, gadsflagmap=GADS_FLAG_MAP,
desired_flags=DESIRED_FLAGS):
desired_flags=DESIRED_FLAGS, extend_reset_flag=False):
"""
Internal function that takes the flags coming from the fits files and turns them into the
bitmasks needed for the HDF file"""
Expand All @@ -505,6 +508,11 @@ def _process_akari_flags(self, pixflag_array, frameflag_array, statusflag_array,
for flag_name, target_bit in desired_flags[flagtype].items():
flag_idx = flagmap[flag_name]
curr_mask = flag_arr[:, flag_idx]
if extend_reset_flag and flagtype == 'pix_flags' and flag_name == 'reset':
# If this option is set, we want to extend the reset flag
# to cover four samples rather than one.
for i in range(1, 4):
curr_mask[i:] = curr_mask[i:] | curr_mask[:-i]
outflags[curr_mask] += 2 ** int(target_bit)
modeflags = desired_flags['mode']
# These numbers are the modes of the pixflag_array field, listed in the AKARI fits data user
Expand Down Expand Up @@ -540,13 +548,13 @@ def inner(self, *args):
@get_chunk
def get_chunk_data(self, band: str, det:str):
# The third entry is just zeroes - i.e. currently we're not operating with a 'psi'.
if self.should_compress:
tod = self.all_chunk_data[band][det]['todz']
else:
tod = self.all_chunk_data[band][det]['tod']
tod = self.all_chunk_data[band][det]['tod']
todz = self.all_chunk_data[band][det]['todz']
psi_arr = np.zeros_like(self.all_chunk_data[band][det]['pix'], dtype=int)
return (tod,
self.all_chunk_data[band][det]['pix'],
todz,
np.array([self.all_chunk_data[band][det]['pix'],
self.all_chunk_data[band][det]['pix_solarcentric']]),
psi_arr,
self.all_chunk_data[band][det]['flag'])

Expand Down Expand Up @@ -578,8 +586,6 @@ def get_chunk_end_earthpos(self):
def get_chunk_satvel(self):
return [0]*3

def get_should_compress_tods(self):
return self.should_compress

def get_gain(self, band: str, detector: str):
return 1
Expand Down
13 changes: 8 additions & 5 deletions commander3/todscripts/AKARI/fits2hdf.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,15 +13,18 @@
#NUM_SEGMENT_PROCESSES = None
# You can run all together, but it's typically faster to start four different
# instances, one for each band
#BANDS = ['065', '090', '140', '160']
BANDS = ['140']
BANDS = ['065', '090', '140', '160']
EXTEND_RESET_FLAG = True

def run_single_band_write(band):
akari_comm_data_adapter = akari_commander_data_adapter.AKARICommanderDataAdapter(
AKARI_FITS_DIR, NSIDE, bands=[band])
AKARI_FITS_DIR, NSIDE, bands=[band],
extend_reset_flag=EXTEND_RESET_FLAG)
comm_todwriter = CommanderHDFWriter(akari_comm_data_adapter)
comm_todwriter.write_hdf_files(OUTPATH, overwrite=True, num_processes=NUM_SEGMENT_PROCESSES,
bands=[band])
comm_todwriter.write_hdf_files(OUTPATH, overwrite=True,
num_processes=NUM_SEGMENT_PROCESSES,
bands=[band],
include_solarcentric_pixpointing=True)

def main():
for band in BANDS:
Expand Down