From 95738b72d543093fde2d23f8dcd8db4a3ef30f26 Mon Sep 17 00:00:00 2001 From: John Henry Beale Date: Mon, 15 Jun 2026 10:10:38 +0200 Subject: [PATCH] BS data now found correctly --- src/clara.py | 518 ++++++++++++++++++++++++++++++++++----------------- 1 file changed, 344 insertions(+), 174 deletions(-) diff --git a/src/clara.py b/src/clara.py index 17ec97b..8ea3769 100755 --- a/src/clara.py +++ b/src/clara.py @@ -12,6 +12,7 @@ import sys import time import zmq import re +import glob from pathlib import Path from loguru import logger @@ -23,7 +24,7 @@ import receive_msg #sys.path.insert(0, os.path.expanduser('/sf/cristallina/applications/mx/jungfraujoch_openapi_python/sf_datafiles')) #from sfdata import SFDataFiles #define log file place: -LOG_FILENAME = time.strftime("/sf/cristallina/data/p23083/res/clara_%Y%m.log") +LOG_FILENAME = time.strftime("/sf/cristallina/data/p23523/res/clara_%Y%m.log") logger.add(LOG_FILENAME, level="INFO", rotation="100MB") @@ -72,54 +73,71 @@ def to_json(obj): """ return json.dumps(obj, default=lambda obj: obj.__dict__, indent=4) - -def wait_for_dataset(h5file, dataset_path, timeout=60, check_interval=0.1): +def wait_for_dataset( + filepath, + dataset_path, + timeout=60, + check_interval=0.1 + ): """ - Wait until an HDF5 dataset exists and contains data. - - Parameters - ---------- - h5file : h5py.File - Open HDF5 file handle - dataset_path : str - Path to dataset inside HDF5 file - timeout : float - Max wait time in seconds - check_interval : float - Time between checks in seconds - + Wait until an HDF5 dataset becomes readable. + + Handles: + - incomplete HDF5 metadata + - delayed writer flushes + - unstable object headers + - live DAQ writing + Returns ------- - dataset - The HDF5 dataset object - - Raises - ------ - TimeoutError - If dataset does not become available in time + np.ndarray + Dataset contents """ - + start = time.time() + attempt = 0 while True: - try: - if dataset_path in h5file: - ds = h5file[dataset_path] - # Optional: ensure dataset is populated - if ds.shape[0] > 0: - return ds + h5file = None + + try: + + h5file = h5py.File(filepath, "r") + + data = h5file[dataset_path][:] + + h5file.close() + + return data except Exception as e: - logger.info(f"Waiting for dataset {dataset_path}: {e}") + + attempt += 1 + + if attempt % 100 == 0: + + logger.info( + f"Still waiting for dataset " + f"{dataset_path} " + f"in {filepath}: {e}" + ) + + try: + h5file.close() + except: + pass if time.time() - start > timeout: + raise TimeoutError( - f"Timeout waiting for dataset {dataset_path}" + f"Timeout waiting for dataset " + f"{dataset_path}" ) time.sleep(check_interval) + # --------class with functions----- @@ -443,202 +461,340 @@ class CollectedH5: return output_file_path + def find_matching_bsdata(self, pulseids_JF, data_dir): + """ + Search all BSDATA files and find pulse ID matches + with the current Jungfrau acquisition. + """ + + matched_results = [] + + bsdata_files = sorted( + glob.glob( + os.path.join(data_dir, "acq*.BSDATA.h5") + ) + ) + + logger.info( + f"Searching {len(bsdata_files)} BSDATA files" + ) + + for bsfile in bsdata_files: + + logger.info(f"Checking BSDATA file: {bsfile}") + + try: + + bsdata = h5py.File(bsfile, "r") + + pulseids_BS = wait_for_dataset( + bsdata, + "/SAR-CVME-TIFALL6:EvtSet/pulse_id", + timeout=60, + check_interval=0.1 + )[:] + + evt_set = wait_for_dataset( + bsdata, + "/SAR-CVME-TIFALL6:EvtSet/data", + timeout=60, + check_interval=0.1 + )[:] + + inters, ind_JF, ind_BS = np.intersect1d( + pulseids_JF, + pulseids_BS, + return_indices=True + ) + + logger.info( + f"{os.path.basename(bsfile)} " + f"matched {len(inters)} pulse IDs" + ) + + if len(inters) > 0: + + matched_results.append( + ( + evt_set, + ind_JF, + ind_BS, + bsfile + ) + ) + + bsdata.close() + + except Exception as e: + + logger.info( + f"Could not process BSDATA file " + f"{bsfile}: {e}" + ) + + logger.info( + f"Found {len(matched_results)} " + f"matching BSDATA files" + ) + + return matched_results + def create_list_file(self): """ Function to generate a list file with the path of the input H5 file :return:None """ + + import glob + print('writing list files') - #Assemble path for raw data + + # Assemble path for raw data pgroup = str(self.message["experiment_group"]) raw = "raw" data = "data" filen = str(self.message["filename"]) + file_path = slash / filen # write to cell file in output folder run_number = self.message["run_number"] acq_number = self.message["file_number"] + name = f'run{run_number:04}_acq{acq_number:04}' - #f = open(name + ".list", "w") - #f.write(str(file_path)) - #f.close() - bsdata_name = Path(filen[:-15]+'.BSDATA.h5') + jf_path = file_path - bsdata_path = slash / bsdata_name - - logger.info("Waiting for bs_data to exist") - start_t = time.time() - - while not bsdata_path.exists(): - current_t = time.time() - elapsed_t = current_t-start_t - - if elapsed_t > 600: - print("Wait for BSdata exceeds 10 minutes. Exiting loop") - break - - time.sleep(1) - - try: - time.sleep(10) - #bsdata = SFDataFiles(bsdata_path) - bsdata = h5py.File(bsdata_path, "r") #r"/sf/cristallina/data/p21630/raw/run0065-lov_movestop_normal_1/data/acq0001.BSDATA.h5" - except Exception as e: - print(f"didn't open bsdata due to error {e}") #_logger.error(f"Cannot open {data_file} (due to {e})") - return - - pulseids_BS = wait_for_dataset( - bsdata, - "/SAR-CVME-TIFALL6:EvtSet/pulse_id", - timeout=600, - check_interval=0.1 - )[:] - - evt_set = wait_for_dataset( - bsdata, - "/SAR-CVME-TIFALL6:EvtSet/data", - timeout=600, - check_interval=0.1 - )[:] - - jf_path= file_path - master_path = slash / Path(filen[:-19]+'_master.h5') + master_path = slash / Path(filen[:-19] + '_master.h5') print(master_path) file_exists = os.path.exists(master_path) st = time.perf_counter() while not file_exists: - time.sleep(1) - file_exists = os.path.exists(master_path) - if time.perf_counter()-st > 1800: - print('Wait for acq master exceeded 30 minutes!!!') - break + time.sleep(1) + file_exists = os.path.exists(master_path) + if time.perf_counter() - st > 1800: + print('Wait for acq master exceeded 30 minutes!!!') + break + + try: - try: print(f"JF file path = {jf_path}, file_path = {file_path}") + print(os.path.exists(jf_path)) - x = h5py.File(jf_path, "r") - #data = SFDataFiles(bsdata_path, jf_path) + + pulseids_JF = wait_for_dataset( + jf_path, + "/entry/xfel/pulseID", + timeout=600, + check_interval=0.1 + ) + except Exception as e: - print(f"didn't open {jf_path} due to error {e}") #_logger.error(f"Cannot open {data_file} (due to {e})") + + print(f"didn't open {jf_path} due to error {e}") + return - pulseids_JF = wait_for_dataset( - x, - "/entry/xfel/pulseID", - timeout=600, - check_interval=0.1 - )[:] - - n_pulse_id = len(pulseids_JF) #- maybe not needed ? - index_dark = [] - index_light = [] - blanks = [] - #print(evt_set) + index_dark = set() + index_light = set() + blanks = set() - _inters, ind_JF, ind_BS = np.intersect1d(pulseids_JF, pulseids_BS, return_indices=True) - #evt_set = evt_set[ind_BS] - #print(len(evt_set)) - - #print(min(ind_BS),max(ind_BS),len(ind_BS)) + # ----------------------------------------------------------------- + # Wait for expected BSDATA file for this acquisition + # ----------------------------------------------------------------- - for idx_JF, idx_BS in zip(ind_JF, ind_BS): - events = evt_set[idx_BS] - - trigger_event = int(self.message["user_data"]["trigger_event"]) + data_dir = os.path.dirname(str(file_path)) - if trigger_event == 63: - - if events[trigger_event] and events[200]: - index_dark.append(idx_JF) + bsdata_name = Path( + filen[:-15] + '.BSDATA.h5' + ) - elif events[200]: - index_light.append(idx_JF) + bsdata_path = slash / bsdata_name - else: - blanks.append(idx_JF) + logger.info( + f"Waiting for BSDATA file: {bsdata_path}" + ) - else: + start_t = time.time() - if events[trigger_event] and events[200]: - index_light.append(idx_JF) + while not bsdata_path.exists(): - elif events[200]: - index_dark.append(idx_JF) + elapsed_t = time.time() - start_t - else: - blanks.append(idx_JF) + if elapsed_t > 600: + logger.info( + "Wait for BSDATA exceeded 10 minutes" + ) - #for i in range(n_pulse_id): + break - # p = pulseids_JF[i] - # q = pulseids_BS[i] - # if p != q: - # #event_i = np.where(pulseids_JF == q)[0] for normal detector - # event_i = np.where(pulseids_BS == p)[0][0]#[::2] #[0][0] - # #event_i = i - # else: - # event_i=i# + time.sleep(0.1) - # events=evt_set[event_i] - #print(events) - #print(len(events)) - #print(event_i) - #print(p,q) - #print(i) + logger.info( + "Expected BSDATA file appeared" + ) - # trigger_event = int(self.message["user_data"]["trigger_event"]) + # ----------------------------------------------------------------- + # Search ALL BSDATA files for matching pulse IDs + # ----------------------------------------------------------------- - # if self.message["user_data"]["trigger_flag"]: - # if events[trigger_event] and events[200]: - # index_light.append(i) + bsdata_files = sorted( + glob.glob( + os.path.join(data_dir, "acq*.BSDATA.h5") + ) + ) - # elif events[200]: - # index_dark.append(i) + logger.info( + f"Found {len(bsdata_files)} BSDATA files" + ) - # else: - # blanks.append(i) + for bsfile in bsdata_files: + + logger.info(f"Checking BSDATA file: {bsfile}") + + try: + + pulseids_BS = wait_for_dataset( + bsfile, + "/SAR-CVME-TIFALL6:EvtSet/pulse_id", + timeout=30, + check_interval=0.1 + )[:] + + evt_set = wait_for_dataset( + bsfile, + "/SAR-CVME-TIFALL6:EvtSet/data", + timeout=30, + check_interval=0.1 + )[:] + + # --------------------------------------------------------- + # Match pulse IDs + # --------------------------------------------------------- + + inters, ind_JF, ind_BS = np.intersect1d( + pulseids_JF, + pulseids_BS, + return_indices=True + ) + + logger.info( + f"{os.path.basename(bsfile)} " + f"matched {len(inters)} pulse IDs" + ) + + trigger_event = int( + self.message["user_data"]["trigger_event"] + ) + + # --------------------------------------------------------- + # Sort light / dark / blank + # --------------------------------------------------------- + + for idx_JF, idx_BS in zip(ind_JF, ind_BS): + + events = evt_set[idx_BS] + + if trigger_event == 63: + + if events[trigger_event] and events[200]: + + index_dark.add(idx_JF) + + elif events[200]: + + index_light.add(idx_JF) + + else: + + blanks.add(idx_JF) + + else: + + if events[trigger_event] and events[200]: + + index_light.add(idx_JF) + + elif events[200]: + + index_dark.add(idx_JF) + + else: + + blanks.add(idx_JF) + + except Exception as e: + + logger.info( + f"Could not process BSDATA file " + f"{bsfile}: {e}" + ) - bsdata.close() - x.close() - acq_dark = [] acq_light = [] - acq_blank = [] + acq_blank = [] + delim = "//" - + if index_dark: - for frame_number in index_dark: - acq_dark.append(f"{jf_path} {delim}{frame_number}") - - file_off = name+'_off.list' - - logger.info(f"List of dark frames : {file_off} , {len(index_dark)} frames") - + + for frame_number in sorted(index_dark): + + acq_dark.append( + f"{jf_path} {delim}{frame_number}" + ) + + file_off = name + '_off.list' + + logger.info( + f"List of dark frames : " + f"{file_off} , {len(index_dark)} frames" + ) + with open(file_off, "w") as f_list: + for frame in acq_dark: - print(f"{frame}", file = f_list) + + print(f"{frame}", file=f_list) if index_light: - for frame_number in index_light: - acq_light.append(f"{jf_path} {delim}{frame_number}") - file_on = name+'_on.list' - logger.info(f"List of light frames : {file_on} , {len(index_light)} frames") + + for frame_number in sorted(index_light): + + acq_light.append( + f"{jf_path} {delim}{frame_number}" + ) + + file_on = name + '_on.list' + + logger.info( + f"List of light frames : " + f"{file_on} , {len(index_light)} frames" + ) + with open(file_on, "w") as f_list: + for frame in acq_light: - print(f"{frame}", file = f_list) + + print(f"{frame}", file=f_list) if blanks: - for frame_number in blanks: - acq_blank.append(f"{jf_path} {delim}{frame_number}") - file_blank = name+'_blank.list' - with open(file_blank, "w") as f_list: - for frame in acq_blank: - print(f"{frame}", file = f_list) - + for frame_number in sorted(blanks): + + acq_blank.append( + f"{jf_path} {delim}{frame_number}" + ) + + file_blank = name + '_blank.list' + + with open(file_blank, "w") as f_list: + + for frame in acq_blank: + + print(f"{frame}", file=f_list) + return None def create_slurm_script(self,trigger, geometry_file_path): @@ -766,12 +922,25 @@ class CollectedH5: # sub.run(["sbatch", "run_SLURM"]) try: - slurm_out = sub.run(["sbatch", "run_SLURM_" + trigger ], capture_output=True) - time.sleep(1) - txt = slurm_out.stdout.decode().split() - # grep the slurm number + slurm_out = sub.run( + ["sbatch", "run_SLURM_" + trigger], + capture_output=True, + text=True + ) + + logger.info(f"sbatch stdout: {slurm_out.stdout}") + logger.info(f"sbatch stderr: {slurm_out.stderr}") + + txt = slurm_out.stdout.split() + + if len(txt) == 0: + raise RuntimeError( + f"sbatch returned empty stdout. stderr={slurm_out.stderr}" + ) + logger.info(f"submitted batch job number: {txt[-1]}") self.message["SlurmJobID_" + trigger] = str(txt[-1]) + except Exception as e: logger.info("Could not submit SLURM job: {}".format(e)) @@ -855,7 +1024,8 @@ if __name__ == "__main__": context = zmq.Context() subscriber = context.socket(zmq.SUB) #subscriber.connect("tcp://sf-broker-01.psi.ch:5555") - subscriber.connect("tcp://sf-daqtest-01.psi.ch:5401") + #subscriber.connect("tcp://sf-daqtest-01.psi.ch:5401") + subscriber.connect("tcp://sf-daq-2.psi.ch:5401") subscriber.setsockopt_string(zmq.SUBSCRIBE, "")