scdms_soudan_parser.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from construct import *\n",
    "\n",
    "format_word = Struct(\n",
    "    \"daq_major\" / Byte,\n",
    "    \"daq_minor\" / Byte,\n",
    "    \"data_format_major\" / Byte,\n",
    "    \"data_format_minor\" / Byte\n",
    ")\n",
    "\n",
    "two_word_file_header = Struct(\n",
    "    \"endian_indicator\" / Int32ul,\n",
    "    \"data_format\" / format_word\n",
    ")\n",
    "\n",
    "detector_hdr = Struct(\n",
    "    \"header_number\" / Int32ul,\n",
    "    \"config_record_len\" / Int32ul,\n",
    "    \"repeat_value\" / Computed(\n",
    "        lambda this: (this.config_record_len // 72) + (this.config_record_len // 144)\n",
    "    )\n",
    ")\n",
    "\n",
    "charge_config_header = Struct(\n",
    "    \"charge_config_len\" / Int32ul,\n",
    "    \"detector_code\" / Int32sl,\n",
    "    \"tower_number\" / Int32sl,\n",
    "    \"channel_post_amp\" / Int32sl,\n",
    "    \"channel_bias\" / Int32sl,\n",
    "    \"rtf_offset\" / Int32sl,\n",
    "    \"delta_t\" / Int32sl,\n",
    "    \"trigger_time\" / Int32sl,\n",
    "    \"trace_len\" / Int32sl\n",
    ")\n",
    "\n",
    "phonon_config_header = Struct(\n",
    "    \"phonon_config_len\" / Int32ul,\n",
    "    \"detector_code\" / Int32sl,\n",
    "    \"tower_number\" / Int32sl,\n",
    "    \"post_amp_gain\" / Int32sl,\n",
    "    \"qet_bias\" / Int32sl,\n",
    "    \"squid_bias\" / Int32sl,\n",
    "    \"squid_lockpoint\" / Int32sl,\n",
    "    \"rtf_offset\" / Int32sl,\n",
    "    \"variable_gain\" / Int32sl,\n",
    "    \"delta_t\" / Int32sl,\n",
    "    \"trigger_time\" / Int32sl,\n",
    "    \"trace_len\" / Int32sl\n",
    ")\n",
    "\n",
    "header_list = Struct(\n",
    "    \"header_number\" / Int32ul,\n",
    "    \"charge_config\" / If(\n",
    "        lambda this: this.header_number == 0x10002,\n",
    "        charge_config_header\n",
    "    ),\n",
    "    \"phonon_config\" / If(\n",
    "        lambda this: this.header_number == 0x10001,\n",
    "        phonon_config_header\n",
    "    )\n",
    ")\n",
    "\n",
    "event_header = Struct(\n",
    "    \"event_header_word\" / Int32ul,\n",
    "    \"event_size\" / Int32ul,\n",
    "    \"event_identifier\" / Computed(\n",
    "        lambda this: (this.event_header_word >> 16) & 0xFFFF\n",
    "    ),\n",
    "    # 0x0: Raw, 0x1: Processed, 0x2: Monte Carlo\n",
    "    \"event_class\" / Computed(\n",
    "        lambda this: (this.event_header_word >> 8) & 0xF\n",
    "    ),\n",
    "    # 0x0: Per Trigger, 0x1: Occasional, 0x2: Begin File Series, 0x3: Begin File\n",
    "    # 0x4: End File, 0x5: End File Series, 0x6: Per Trigger w/ Detectors that Cross Threshold\n",
    "    \"event_category\" / Computed(\n",
    "        lambda this: (this.event_header_word >> 12) & 0xF\n",
    "    ),\n",
    "    # 0x0: Wimp Search, 0x1: 60Co Calibration, 0x2: 60Co Low Energy Calibration,\n",
    "    # 0x3: Neutron Calibration, 0x4: Random Triggers, 0x5: Pulse Triggers\n",
    "    # 0x6: Test, 0x7: Data Monitering Event, 0x8: 137Cs Calibration\n",
    "    \"event_type\" / Computed(\n",
    "        lambda this: (this.event_header_word & 0xFF)\n",
    "    )\n",
    ")\n",
    "\n",
    "administrative_record = Struct(\n",
    "    \"admin_header\" / Int32ul,\n",
    "    \"admin_len\" / Int32ul,\n",
    "    \"series_number_1\" / Int32ul,\n",
    "    \"series_number_2\" / Int32ul,\n",
    "    \"event_number_in_series\" / Int32ul,\n",
    "    \"seconds_from_epoch\" / Int32ul,\n",
    "    # Epoch defined as Jan 1st 1904 for SUF (MAC Artifact)\n",
    "    # Epoch defined as Jan 1st 1970 for Soudan\n",
    "    \"time_from_last_event\" / Int32ul,\n",
    "    \"live_time_from_last_event\" / Int32ul\n",
    ")\n",
    "\n",
    "trace_record = Struct(\n",
    "    \"trace_header\" / Int32ul,\n",
    "    \"trace_len\" / Int32ul,\n",
    "    \"trace_bookkeeping_header\" / Int32ul,\n",
    "    \"bookkeeping_len\" / Int32ul,\n",
    "    \"digitizer_base_address\" / Int32ul,\n",
    "    \"digitizer_channel\" / Int32ul,\n",
    "    \"detector_code\" / Int32ul,\n",
    "    \"timebase_header\" / Int32ul,\n",
    "    \"timebase_len\" / Int32ul,\n",
    "    \"t0_in_ns\" / Int32ul,\n",
    "    \"delta_t_ns\" / Int32ul,\n",
    "    \"num_of_points\" / Int32ul,\n",
    "    \"second_trace_header\" / Int32ul,\n",
    "    \"num_samples\" / Int32ul\n",
    "    # Should be a power of two (1024, 2048, etc)\n",
    ")\n",
    "\n",
    "data_sample = Struct(\n",
    "    \"data_selection\" / Int32ul,\n",
    "    \"sample_a\" / Computed(\n",
    "        lambda this: (this.data_selection >> 16) & 0xFFFF\n",
    "    ),\n",
    "    \"sample_b\" / Computed(\n",
    "        lambda this: (this.data_selection & 0xFFFF)\n",
    "    )\n",
    ")\n",
    "\n",
    "\n",
    "trace_data = Struct(\n",
    "    \"trace_rcrds\" / trace_record,\n",
    "    \"sample_data\" / Array(\n",
    "        this.trace_rcrds.num_samples // 2,\n",
    "        data_sample\n",
    "    )\n",
    ")\n",
    "\n",
    "soudan_history_buffer = Struct(\n",
    "    \"history_buffer_header\" / Int32ul,\n",
    "    \"history_buffer_len\" / Int32ul,\n",
    "    \"num_time_nvt\" / Int32ul,\n",
    "    \"time_nvt\" / Array(\n",
    "        this.num_time_nvt,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"num_veto_mask_words\" / Int32ul,\n",
    "    \"time_n_minus_veto_mask\" / Array(\n",
    "        this.num_time_nvt * this.num_veto_mask_words,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"num_trigger_times\" / Int32ul,\n",
    "    \"trigger_times\" / Array(\n",
    "        this.num_trigger_times,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"num_trigger_mask_words\" / Int32ul,\n",
    "    \"trig_times_minus_trig_mask\" / Array(\n",
    "        this.num_trigger_times * this.num_trigger_mask_words,\n",
    "        Int32ul\n",
    "    )\n",
    ")\n",
    "\n",
    "trigger_record = Struct(\n",
    "    \"trigger_header\" / Int32ul,\n",
    "    \"trigger_len\" / Int32ul,\n",
    "    \"trigger_time\" / Int32ul,\n",
    "    \"individual_trigger_masks\" / Array(\n",
    "        6,\n",
    "        Int32ul\n",
    "    )\n",
    ")\n",
    "\n",
    "tlb_trigger_mask_record = Struct(\n",
    "    \"tlb_mask_header\" / Int32ul,\n",
    "    \"tlb_len\" / Int32ul,\n",
    "    \"tower_mask\" / Array(\n",
    "        6,\n",
    "        Int32ul\n",
    "    )\n",
    ")\n",
    "\n",
    "gps_data = Struct(\n",
    "    \"tlb_mask_header\" / Int32ul,\n",
    "    \"length\" / Int32ul,\n",
    "    \"gps_year_day\" / If(\n",
    "        this.length > 0,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"gps_status_hour_minute_second\" / If(\n",
    "        this.length > 0,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"gps_microsecs_from_gps_second\" / If(\n",
    "        this.length > 0,\n",
    "        Int32ul\n",
    "    )\n",
    ")\n",
    "\n",
    "detector_trigger_threshold_data = Struct(\n",
    "    \"threshold_header\" / Int32ul,\n",
    "    \"len_to_next_header\" / Int32ul,\n",
    "    \"minimum_voltage_level\" / Int32ul,\n",
    "    \"maximum_voltage_level\" / Int32ul,\n",
    "    \"dynamic_range\" / Int32ul,\n",
    "    \"tower_number\" / Int32ul,\n",
    "    \"detector_codes\" / Array(\n",
    "        6,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"operations_codes\" / Array(\n",
    "        9,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"adc_values\" / Array(\n",
    "        54,\n",
    "        Int32ul\n",
    "    )\n",
    ")\n",
    "\n",
    "detector_trigger_rates = Struct(\n",
    "    \"detector_trigger_header\" / Int32ul,\n",
    "    \"len_to_next_header\" / Int32ul,\n",
    "    \"clocking_interval\" / Int32ul,\n",
    "    \"tower_number\" / Int32ul,\n",
    "    \"detector_codes\" / Array(\n",
    "        6,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"j_codes\" / Array(\n",
    "        5,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"counter_values\" / Array(\n",
    "        30,\n",
    "        Int32ul\n",
    "    )\n",
    ")\n",
    "\n",
    "veto_trigger_rates = Struct(\n",
    "    \"veto_trigger_header\" / Int32ul,\n",
    "    \"len_to_next_header\" / Int32ul,\n",
    "    \"clocking_interval\" / Int32ul,\n",
    "    \"num_entries\" / Int32ul,\n",
    "    \"detector_code\" / Array(\n",
    "        this.num_entries,\n",
    "        Int32ul\n",
    "    ),\n",
    "    \"counter_value_det_code\" / Array(\n",
    "        this.num_entries,\n",
    "        Int32ul\n",
    "    )\n",
    ")\n",
    "\n",
    "\n",
    "logical_records = Struct(\n",
    "    \"event_hdr\" / Peek(Int32ul),  # Peek to check first\n",
    "    \"next_section\" / Struct(\n",
    "        \"next_header\" / Peek(Int32ul),  # Peek without consuming\n",
    "        \"section\" / Switch(\n",
    "            lambda this: (\n",
    "                this.next_header if ((this.next_header >> 16) != 0xA980) \n",
    "                else 0xA980  # Use 0xA980 as identifier for event_header\n",
    "            ),\n",
    "            {\n",
    "                0xA980: event_header,\n",
    "                0x00000002: administrative_record,\n",
    "                0x00000011: trace_data,\n",
    "                0x00000021: soudan_history_buffer,\n",
    "                0x00000060: gps_data,\n",
    "                0x00000080: trigger_record,\n",
    "                0x00000081: tlb_trigger_mask_record,\n",
    "                0x00000022: detector_trigger_rates,\n",
    "                0x00000031: veto_trigger_rates,\n",
    "            }\n",
    "        )\n",
    "    )\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import h5py\n",
    "\n",
    "def get_detector_code_info(detector_code):\n",
    "    \"\"\"\n",
    "    Detector codes come in the form xxxyyyzzz\n",
    "\n",
    "    xxx - detector type (0-999)\n",
    "\n",
    "    yyy - detector number (0-999)\n",
    "\n",
    "    zzz - channel number (0-999)\n",
    "\n",
    "    Detector types and channel numbers should not exceed two digits.\n",
    "    \"\"\"\n",
    "    detector_type_dictionary = [\n",
    "        # Channel numbers gives charge or phonon\n",
    "        {'ID': 1,  'Det Type': 'Blip', 'Charge': [1,2], 'Phonon': [3,4]},\n",
    "        {'ID': 2,  'Det Type': 'Flip', 'Charge': [0,1], 'Phonon': [2,3,4,5]},\n",
    "        {'ID': 3,  'Det Type': 'Veto', 'Charge': 'None', 'Phonon': 'None'},\n",
    "        {'ID': 4,  'Det Type': 'ZIP', 'Charge': [0,1], 'Phonon': [2,3,4,5]},\n",
    "        {'ID': 5,  'Det Type': 'mercedes ZIP', 'Charge': [0,1], 'Phonon': [2,3,4,5]},\n",
    "        {'ID': 6,  'Det Type': 'endcap (class I)', 'Charge': [0,1], 'Phonon': [2,3,4,5]},\n",
    "        {'ID': 7,  'Det Type': 'endcap (class II)', 'Charge': [0], 'Phonon': [1,2]},\n",
    "        {'ID': 10, 'Det Type': 'iZIP (class I)', 'Charge': [0,1,6,7], 'Phonon': [2,3,4,5,8,9,10,11]},\n",
    "        {'ID': 11, 'Det Type': 'iZIP (class II)', 'Charge': [0,1,6,7], 'Phonon': [2,3,4,5,8,9,10,11]}\n",
    "    ]\n",
    "\n",
    "    if len(str(detector_code)) == 8:\n",
    "        detector_type = int(str(detector_code)[:2])\n",
    "    else:\n",
    "        detector_type = int(str(detector_code)[:1])\n",
    "    detector_number = int(str(detector_code)[-6:-3]) # Not used here, returned for utility\n",
    "    channel_number = int(str(detector_code)[-3:])\n",
    "\n",
    "    try:\n",
    "        detector_info = next((d for d in detector_type_dictionary if d['ID'] == detector_type))\n",
    "    except:\n",
    "        #print(f'Error on detector code: {detector_code}\\nDetector type not found.')\n",
    "        det_type = 'Error'\n",
    "        charge = False\n",
    "        phonon = False\n",
    "        veto = False\n",
    "        error = True\n",
    "    \n",
    "    try:\n",
    "        det_type = detector_info['Det Type']\n",
    "        # Reverse lookup for the channel number\n",
    "        if isinstance(detector_info['Charge'], list) and channel_number in detector_info['Charge']:\n",
    "            charge = True\n",
    "            phonon = False\n",
    "            veto = False\n",
    "            error = False\n",
    "        elif isinstance(detector_info['Phonon'], list) and channel_number in detector_info['Phonon']:\n",
    "            charge = False\n",
    "            phonon = True\n",
    "            veto = False\n",
    "            error = False\n",
    "        elif detector_info['Det Type'] == 'Veto':\n",
    "            charge = False\n",
    "            phonon = False\n",
    "            veto = True\n",
    "            error = False\n",
    "        else:\n",
    "            charge = False\n",
    "            phonon = False\n",
    "            veto = False\n",
    "            error = True\n",
    "    except:\n",
    "        det_type = 'Error'\n",
    "        charge = False\n",
    "        phonon = False\n",
    "        veto = False\n",
    "        error = True\n",
    "\n",
    "    return det_type, charge, phonon, veto, error, detector_number\n",
    "\n",
    "# Test detector code\n",
    "type, charge, phonon, veto, error, detector_number = get_detector_code_info(3555001)\n",
    "print(f'Type: {type}, Charge: {charge}, Phonon: {phonon}, Veto: {veto}, Error: {error}, Detector Number: {detector_number}')\n",
    "\n",
    "# Print structure function to inspect groups mid processing\n",
    "def print_structure(name, obj):\n",
    "    if isinstance(obj, h5py.Group):\n",
    "        print(f'Group: {name}')\n",
    "    elif isinstance(obj, h5py.Dataset):\n",
    "        print(f\"Dataset: {name}, Shape: {obj.shape}, Data type: {obj.dtype}\")\n",
    "        # Print the dataset values (for small datasets)\n",
    "        data = obj[()]\n",
    "        print(f\"Values: {data}\")\n",
    "\n",
    "soudan = Struct(\n",
    "    \"file_hdr\" / two_word_file_header,\n",
    "    \"detector_hdr\" / detector_hdr,\n",
    "    \"hdrs\" / Array(\n",
    "        this._root.detector_hdr.repeat_value,\n",
    "        header_list\n",
    "    ),\n",
    "    \"logical_rcrds\" / GreedyRange(logical_records)\n",
    ")\n",
    "\n",
    "# Parsing instructions to create smaller files faster\n",
    "test = Struct(\n",
    "    \"file_hdr\" / two_word_file_header,\n",
    "    \"detector_hdr\" / detector_hdr,\n",
    "    \"hdrs\" / Array(\n",
    "        this._root.detector_hdr.repeat_value,\n",
    "        header_list\n",
    "    ),\n",
    "    \"logical_rcrds\" / Array(\n",
    "        1500,\n",
    "        logical_records\n",
    "    )\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open(\"minimal_sample.soudan\", 'rb') as input_f:\n",
    "    raw_data = input_f.read()\n",
    "    # Switch between soudan and test for different amounts of parsing\n",
    "    parsed_data = soudan.parse(raw_data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "parsed_data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def parse_file(input_path, output_path, use_test_parse=True):\n",
    "    if use_test_parse:\n",
    "        print(f'Parsing {input_path}...')\n",
    "    with open(input_path, 'rb') as input_f:\n",
    "        raw_data = input_f.read()\n",
    "        # Switch between soudan and test for different amounts of parsing\n",
    "        if use_test_parse:\n",
    "            parsed_data = test.parse(raw_data)\n",
    "        else:\n",
    "            parsed_data = soudan.parse(raw_data)\n",
    "\n",
    "    with h5py.File(output_path, 'w') as output_f:\n",
    "        \n",
    "        # Initializing header groups to fill with datasets\n",
    "        file_hdr_grp = output_f.create_group('file_hdr')\n",
    "        detector_hdr_grp = output_f.create_group('detector_hdr')\n",
    "\n",
    "        # Initializing arrays for the header information\n",
    "        file_hdr_word_list = []\n",
    "\n",
    "        # file_hdr and detector_hdr contain no arrays\n",
    "        if use_test_parse:\n",
    "            print('Parsing file headers...')\n",
    "        for file_hdr_type in parsed_data.file_hdr:\n",
    "            hdr_type_grp = file_hdr_grp.create_group(f'{file_hdr_type}')\n",
    "            file_hdr_word_list.append(hdr_type_grp)\n",
    "            if file_hdr_type == 'data_format':\n",
    "                for attr_name in ['daq_major', 'daq_minor', 'data_format_major', 'data_format_minor']:\n",
    "                    if hasattr(parsed_data.file_hdr.data_format, attr_name):\n",
    "                        attr_value = getattr(parsed_data.file_hdr.data_format, attr_name)\n",
    "                        hdr_type_grp.create_dataset(attr_name, data=attr_value)\n",
    "            elif file_hdr_type == \"endian_indicator\":\n",
    "                hdr_type_grp.create_dataset('endian_indicator', data=parsed_data.file_hdr.endian_indicator)\n",
    "        if use_test_parse:\n",
    "            print(f'Done.\\nParsing detector headers...')\n",
    "        for det_data_type in parsed_data.detector_hdr:\n",
    "            if det_data_type == 'header_number':\n",
    "                detector_hdr_grp.create_dataset('header_number', data=parsed_data.detector_hdr.header_number)\n",
    "            elif det_data_type == 'config_record_len':\n",
    "                detector_hdr_grp.create_dataset('config_record_len', data=parsed_data.detector_hdr.config_record_len)\n",
    "            elif det_data_type == 'repeat_value':\n",
    "                detector_hdr_grp.create_dataset('repeat_value', data=parsed_data.detector_hdr.repeat_value)\n",
    "\n",
    "\n",
    "        # hdrs contains an array of charge and phonon headers\n",
    "        hdrs_grp           = output_f.create_group('hdrs')\n",
    "        charge_config_grp  = hdrs_grp.create_group('charge_config')\n",
    "        phonon_config_grp  = hdrs_grp.create_group('phonon_config')\n",
    "        hdrs_array         = []\n",
    "        charge_config_list = []\n",
    "        phonon_config_list = []\n",
    "\n",
    "        charge_hdr_count = 0\n",
    "        phonon_hdr_count = 0\n",
    "\n",
    "        # Create groups for each header and populate them with relevant datasets\n",
    "        if use_test_parse:\n",
    "            print(f'Done.\\nCreating charge and phonon header groups...')\n",
    "        for i, header in enumerate(parsed_data.hdrs):\n",
    "            # Collecting charge_config data\n",
    "            if header.header_number == 0x10002:\n",
    "                # HDF5 groups require unique names if at same level of structure\n",
    "                charge_config_hdr_grp = charge_config_grp.create_group(f'charge_config_{charge_hdr_count}')\n",
    "                charge_hdr_count += 1\n",
    "                charge_config_list.append(charge_config_hdr_grp)\n",
    "                hdrs_array.append(charge_config_hdr_grp)\n",
    "                for attr_name in ['charge_config_len', 'detector_code', 'tower_number',\n",
    "                                  'channel_post_amp', 'rtf_offset', 'delta_t', 'trigger_time',\n",
    "                                  'trace_len']:\n",
    "                    if hasattr(header.charge_config, attr_name):\n",
    "                        attr_value = getattr(header.charge_config, attr_name)\n",
    "                        charge_config_hdr_grp.create_dataset(attr_name, data=attr_value)\n",
    "                \n",
    "            # Collecting phonon_config data\n",
    "            elif header.header_number == 0x10001:\n",
    "                phonon_config_hdr_grp = phonon_config_grp.create_group(f'phonon_config_{phonon_hdr_count}')\n",
    "                phonon_hdr_count += 1\n",
    "                phonon_config_list.append(header)\n",
    "                hdrs_array.append(phonon_config_hdr_grp)\n",
    "                for attr_name in ['phonon_config_len', 'detector_code', 'tower_number',\n",
    "                                  'post_amp_gain', 'qet_bias', 'squid_bias', 'squid_lockpoint',\n",
    "                                  'rtf_offset', 'variable_gain', 'delta_t', 'trigger_time', 'trace_len']:\n",
    "                    if hasattr(header.phonon_config, attr_name):\n",
    "                        attr_value = getattr(header.phonon_config, attr_name)\n",
    "                        phonon_config_hdr_grp.create_dataset(attr_name, data=attr_value)\n",
    "        \n",
    "        # Dictionary for creating groups\n",
    "        logical_record_options = {\n",
    "            0xA980: \"event_header\",\n",
    "            0x00000002: \"admin_rcrd\",\n",
    "            0x00000011: \"trace_data\",\n",
    "            0x00000021: \"soudan_history_buffer\",\n",
    "            0x00000060: \"gps_data\",\n",
    "            0x00000080: \"trigger_rcrd\",\n",
    "            0x00000081: \"tlb_trigger_mask_rcrd\",\n",
    "            0x00000022: \"detector_trigger_rates\",\n",
    "            0x00000031: \"veto_trigger_rates\",\n",
    "        }\n",
    "\n",
    "        # Creating groups that can hold each event's records\n",
    "        logical_rcrd_grp       = output_f.create_group('logical_rcrds')\n",
    "        event_hdr_grp          = logical_rcrd_grp.create_group('event_hdr')\n",
    "        admin_rcrd_grp         = logical_rcrd_grp.create_group('admin_rcrd')\n",
    "        trigger_rcrd_grp       = logical_rcrd_grp.create_group('trigger_rcrd')\n",
    "        tlb_trig_mask_rcrd_grp = logical_rcrd_grp.create_group('tlb_trig_mask_rcrd')\n",
    "        gps_data_grp           = logical_rcrd_grp.create_group('gps_data')\n",
    "        trace_data_grp         = logical_rcrd_grp.create_group('trace_data')\n",
    "        soudan_buffer_grp      = logical_rcrd_grp.create_group('soudan_buffer')\n",
    "        detector_trig_grp      = logical_rcrd_grp.create_group('detector_trigger_rates')\n",
    "        veto_trig_grp          = logical_rcrd_grp.create_group('veto_trigger_rates')\n",
    "\n",
    "        # Initializing arrays to store the created groups of logical_rcrd data\n",
    "        event_hdr_array          = []\n",
    "        admin_rcrd_array         = []\n",
    "        trigger_rcrd_array       = []\n",
    "        tlb_trig_mask_rcrd_array = []\n",
    "        gps_data_array           = []\n",
    "        trace_data_array         = []\n",
    "        soudan_buffer_array      = []\n",
    "        detector_trig_array      = []\n",
    "        veto_trig_array          = []\n",
    "\n",
    "        # Initializing arrays for storing trace_data data samples\n",
    "        phonon_pulse_array       = []\n",
    "        charge_pulse_array       = []\n",
    "        veto_array               = []\n",
    "        error_array              = []\n",
    "\n",
    "        # Counters for group naming iteration\n",
    "        event_count    = 0\n",
    "        admin_count    = 0\n",
    "        trigger_count  = 0\n",
    "        tlb_count      = 0\n",
    "        gps_count      = 0\n",
    "        trace_count    = 0\n",
    "        soudan_count   = 0\n",
    "        detector_count = 0\n",
    "        veto_count     = 0\n",
    "\n",
    "        # Using to find the unique information related to these groups\n",
    "        trace_set = set()\n",
    "        detector_type_set = set()\n",
    "        detector_code_set = set()\n",
    "        detector_number_set = set()\n",
    "        event_set = set()\n",
    "        charge_set = set()\n",
    "        phonon_set = set()\n",
    "        error_set = set()\n",
    "        veto_set = set()\n",
    "\n",
    "        # Dictionary to hold individual counters for each trace group\n",
    "        group_counters = {}\n",
    "\n",
    "        # Create a group with the series number to store trace and event data\n",
    "        if use_test_parse:\n",
    "            print(f'Done.\\nFinding series number...')\n",
    "        for record_option in parsed_data.logical_rcrds:\n",
    "            for value, type in logical_record_options.items():\n",
    "                if record_option.next_section.next_header == value:\n",
    "                    if type == 'admin_rcrd':\n",
    "                        series_number_1 = record_option.next_section.section.series_number_1\n",
    "                        series_number_2 = record_option.next_section.section.series_number_2\n",
    "                        series_number = f'{series_number_1}{series_number_2}'\n",
    "        series_grp = output_f.create_group(f'S{series_number}')\n",
    "        if use_test_parse:\n",
    "            print(f'Series number: {series_number}')\n",
    "\n",
    "        # Loop through all of the Structs in logical_rcrds\n",
    "        if use_test_parse:\n",
    "            print(f'Parsing logical records...')\n",
    "        for i, record_option in enumerate(parsed_data.logical_rcrds):\n",
    "            # Handle event headers separately\n",
    "            if (record_option.next_section.next_header >> 16) == 0xA980:\n",
    "                # Storing event_hdr data\n",
    "                events = []\n",
    "                # Loop through attributes of event_hdr and store them in event_hdr_grp_i\n",
    "                event_hdr_grp_i = event_hdr_grp.create_group(f'event_group_{event_count}')\n",
    "                event_count += 1\n",
    "                # Skip event_header_word and event_identifier\n",
    "                # Identifier is identical for each event, and header_word\n",
    "                # is split for other information\n",
    "                if use_test_parse:\n",
    "                    print('Event header data:')\n",
    "                for attr_name in ['event_size', 'event_class', 'event_category', 'event_type']:\n",
    "                    if hasattr(record_option.next_section.section, attr_name):\n",
    "                        attr_value = getattr(record_option.next_section.section, attr_name)\n",
    "                        event_data = event_hdr_grp_i.create_dataset(attr_name, data=attr_value)\n",
    "                        events.append(event_data)\n",
    "                        if use_test_parse:\n",
    "                            print(f'{attr_name}: {attr_value}')\n",
    "                event_hdr_array.append(events)\n",
    "            \n",
    "            # Iterate through the options in logical_rcrds using the type dictionary\n",
    "            for value, type in logical_record_options.items():\n",
    "                if record_option.next_section.next_header == value:\n",
    "                    if type == 'admin_rcrd':\n",
    "                        # Store admin_rcrd data in an array\n",
    "                        admins = []\n",
    "                        admin_group_i = admin_rcrd_grp.create_group(f'{type}_group_{admin_count}')\n",
    "                        admin_count += 1\n",
    "                        if use_test_parse:\n",
    "                            print(f'\\nAdministrative record data:')\n",
    "                        for attr_name in ['admin_header', 'admin_len', 'series_number_1', 'series_number_2',\n",
    "                                        'event_number_in_series', 'seconds_from_epoch', 'time_from_last_event',\n",
    "                                        'live_time_from_last_event']:\n",
    "                            if hasattr(record_option.next_section.section, attr_name):\n",
    "                                attr_value = getattr(record_option.next_section.section, attr_name)\n",
    "                                admin_data = admin_group_i.create_dataset(attr_name, data=attr_value)\n",
    "                                admins.append(admin_data)\n",
    "                                if use_test_parse:\n",
    "                                    print(f'{attr_name}: {attr_value}')\n",
    "                        # Store each admin_rcrd data array in higher level array\n",
    "                        admin_rcrd_array.append(admins)\n",
    "\n",
    "                        # Create separate subgroup for series to hold event information\n",
    "                        event_number_i = record_option.next_section.section.event_number_in_series\n",
    "                        event_grp_i = series_grp.create_group(f'E{event_number_i}')\n",
    "                        event_set.add(event_number_i)\n",
    "                \n",
    "                    if type == 'trigger_rcrd':\n",
    "                        # Store trigger_rcrd data in an array\n",
    "                        triggers = []\n",
    "                        trigger_group_i = trigger_rcrd_grp.create_group(f'{type}_group_{trigger_count}')\n",
    "                        trigger_count += 1\n",
    "                        if use_test_parse:\n",
    "                            print(f'\\nTrigger record data:')\n",
    "                        for attr_name in ['trigger_header', 'trigger_len', 'trigger_time', 'individual_trigger_masks']:\n",
    "                            if hasattr(record_option.next_section.section, attr_name):\n",
    "                                attr_value = getattr(record_option.next_section.section, attr_name)\n",
    "                                trigger_data = trigger_group_i.create_dataset(attr_name, data=attr_value)\n",
    "                                triggers.append(trigger_data)\n",
    "                                if use_test_parse:\n",
    "                                    print(f'{attr_name}: {attr_value}')\n",
    "                        # Store each trigger_rcrd array in higher level array\n",
    "                        trigger_rcrd_array.append(triggers)\n",
    "\n",
    "                    if  type == 'tlb_trigger_mask_rcrd':\n",
    "                        # Store tlb_trigger_mask_rcrd data in an array\n",
    "                        tlb_trig_mask = []\n",
    "                        tlb_trig_group_i = tlb_trig_mask_rcrd_grp.create_group(f'{type}_group_{tlb_count}')\n",
    "                        tlb_count += 1\n",
    "                        if use_test_parse:\n",
    "                                    print(f'\\nTlb trigger mask record data:')\n",
    "                        for attr_name in ['tlb_mask_header', 'tlb_len', 'tower_mask']:\n",
    "                            if hasattr(record_option.next_section.section, attr_name):\n",
    "                                attr_value = getattr(record_option.next_section.section, attr_name)\n",
    "                                tlb_trig_data = tlb_trig_group_i.create_dataset(attr_name, data=attr_value)\n",
    "                                tlb_trig_mask.append(tlb_trig_data)\n",
    "                                if use_test_parse:\n",
    "                                    print(f'{attr_name}: {attr_value}')\n",
    "                        # Store each array in higher level array\n",
    "                        tlb_trig_mask_rcrd_array.append(tlb_trig_mask)\n",
    "\n",
    "                    if type == 'gps_data':\n",
    "                        # Storing gps_data in array\n",
    "                        gps = []\n",
    "                        # Loop through attributes of gps_data and store them in gps_data_group_i\n",
    "                        gps_data_group_i = gps_data_grp.create_group(f'{type}_group_{gps_count}')\n",
    "                        gps_count += 1\n",
    "                        if use_test_parse:\n",
    "                                    print(f'\\nGPS data:')\n",
    "                        #gps_data has a variable number of attributes depending on length value\n",
    "                        attributes = ['tlb_mask_header', 'length']\n",
    "                        if record_option.next_section.section.length != 0:\n",
    "                            attributes.extend(['gps_year_day', 'gps_status_hour_minute_second', 'gps_microsecs_from_gps_second'])\n",
    "                        for attr_name in attributes:\n",
    "                            if hasattr(record_option.next_section.section, attr_name):\n",
    "                                attr_value = getattr(record_option.next_section.section, attr_name)\n",
    "                                gps_dataset = gps_data_group_i.create_dataset(attr_name, data=attr_value)\n",
    "                                gps.append(gps_dataset)\n",
    "                                if use_test_parse:\n",
    "                                    print(f'{attr_name}: {attr_value}')\n",
    "                        # Store gps array in higher level array\n",
    "                        gps_data_array.append(gps)\n",
    "\n",
    "                    # trace_data contains the actual data samples taken by the detector\n",
    "                    # As well as the header information describing the events and detectors\n",
    "                    if type == \"trace_data\":\n",
    "                        if use_test_parse:\n",
    "                            print(f'\\nTrace data:')\n",
    "                        if record_option.next_section.section.trace_rcrds:\n",
    "                            trace_record_group_i = trace_data_grp.create_group(f'{type}_group_{trace_count}')\n",
    "                            #trace_record_group_i = event_grp_i.create_group(f'{type}_group_{trace_count}')\n",
    "                            trace_count += 1\n",
    "                            # Store trace_rcrd data in an array\n",
    "                            trace_rcrd = []\n",
    "                            for attr_name in ['trace_header', 'trace_len', 'trace_bookkeeping_header', 'bookkeeping_len',\n",
    "                                  'digitizer_base_address', 'digitizer_channel', 'detector_code', 'timebase_header',\n",
    "                                  'timebase_len', 't0_in_ns', 'delta_t_ns', 'num_of_points', 'second_trace_header',\n",
    "                                  'num_samples']:\n",
    "                                if hasattr(record_option.next_section.section.trace_rcrds, attr_name):\n",
    "                                    attr_value = getattr(record_option.next_section.section.trace_rcrds, attr_name)\n",
    "                                    trace_rcrd_dataset = trace_record_group_i.create_dataset(attr_name, data=attr_value)\n",
    "                                    trace_rcrd.append(trace_rcrd_dataset)\n",
    "                                    if use_test_parse:\n",
    "                                        print(f'{attr_name}: {attr_value}')\n",
    "                            # Store each trace_rcrd array in higher level array\n",
    "                            trace_data_array.append(trace_rcrd)\n",
    "\n",
    "                        if record_option.next_section.section.sample_data:\n",
    "                            # Store data samples in trace array          \n",
    "                            trace = []\n",
    "                            for data in record_option.next_section.section.sample_data:\n",
    "                                trace.append(data.sample_a)\n",
    "                                trace.append(data.sample_b)\n",
    "\n",
    "                        try:\n",
    "                            detector_code = record_option.next_section.section.trace_rcrds.detector_code\n",
    "                            detector_code_set.add(detector_code)\n",
    "                            det_type, charge, phonon, veto, error, detector_number = get_detector_code_info(int(detector_code))\n",
    "                            \n",
    "                            detector_type_set.add(det_type)\n",
    "                            detector_number_set.add(detector_number)\n",
    "                            trace_set.add(len(trace))\n",
    "\n",
    "                            det_group_name = f'det_code_{detector_code}'\n",
    "\n",
    "                            if det_group_name not in group_counters:\n",
    "                                group_counters[det_group_name] = 0\n",
    "\n",
    "                            trace_group_count = group_counters[det_group_name]\n",
    "\n",
    "                            if det_group_name in event_grp_i:\n",
    "                                detector_group = event_grp_i[det_group_name]\n",
    "\n",
    "                                trace_dataset_name = f'trace_{trace_group_count}'\n",
    "                                if trace_dataset_name not in detector_group:\n",
    "                                    detector_group.create_dataset(trace_dataset_name, data=trace)\n",
    "                            \n",
    "                            else:\n",
    "                                try:\n",
    "                                    detector_group = event_grp_i.create_group(f'{det_group_name}')\n",
    "                                    if 'detector_type' not in detector_group:\n",
    "                                        detector_group.create_dataset(f'detector_type', data=f'{det_type}')\n",
    "                                    if 'detector_number' not in detector_group:\n",
    "                                        detector_group.create_dataset(f'detector_number', data=detector_number)\n",
    "\n",
    "                                    if charge:\n",
    "                                        type = 'Charge'\n",
    "                                        charge_pulse_array.append(trace)\n",
    "                                        charge_set.add(len(trace))\n",
    "                                    elif phonon:\n",
    "                                        type = 'Phonon'\n",
    "                                        phonon_pulse_array.append(trace)\n",
    "                                        phonon_set.add(len(trace))\n",
    "                                    elif veto:\n",
    "                                        type = 'Veto'\n",
    "                                        veto_array.append(trace)\n",
    "                                        veto_set.add(len(trace))\n",
    "                                    elif error:\n",
    "                                        type = 'Error'\n",
    "                                        error_array.append(trace)\n",
    "                                        error_set.add(len(trace))\n",
    "                                    \n",
    "                                    detector_group.create_dataset(f'trace', data=trace)\n",
    "                                    group_counters[det_group_name] += 1\n",
    "                                    detector_group.create_dataset(f'trace_type', data=type)\n",
    "                                except Exception as e:\n",
    "                                    print(f'Dataset error:\\n{e}')\n",
    "                        except Exception as e:\n",
    "                            print(f'Detector group error:\\n{e}')\n",
    "\n",
    "                    if type == 'soudan_history_buffer':\n",
    "                        # Store soudan_history_buffer data in an array\n",
    "                        soudan_buffer = []\n",
    "                        soudan_buffer_group_i = soudan_buffer_grp.create_group(f'{type}_group_{soudan_count}')\n",
    "                        soudan_count += 1\n",
    "                        if use_test_parse:\n",
    "                            print(f'\\nSoudan history buffer data:')\n",
    "                        for attr_name in ['history_buffer_header', 'history_buffer_len', 'num_time_nvt', 'time_nvt',\n",
    "                              'num_veto_mask_words', 'time_n_minus_veto_mask', 'num_trigger_times', \n",
    "                              'trigger_times', 'num_trigger_mask_words', 'trig_times_minus_trig_mask']:\n",
    "                            if hasattr(record_option.next_section.section, attr_name):\n",
    "                                attr_value = getattr(record_option.next_section.section, attr_name)\n",
    "                                soudan_buffer_data = soudan_buffer_group_i.create_dataset(attr_name, data=attr_value)\n",
    "                                soudan_buffer.append(soudan_buffer_data)\n",
    "                                if use_test_parse:\n",
    "                                    print(f'{attr_name}: {attr_value}')\n",
    "                        # Store each soudan_history_buffer array in higher level array\n",
    "                        soudan_buffer_array.append(soudan_buffer)\n",
    "\n",
    "                    if type == 'detector_trigger_rates':\n",
    "                        # Store detector_trigger_rate data in an array\n",
    "                        detector = []\n",
    "                        detector_group_i = detector_trig_grp.create_group(f'{type}_group_{detector_count}')\n",
    "                        detector_count += 1\n",
    "                        if use_test_parse:\n",
    "                            print(f'\\nDetector trigger rates data:')\n",
    "                        for attr_name in ['detector_trigger_header', 'len_to_next_header', 'clocking_interval',\n",
    "                                          'tower_number', 'detector_codes', 'j_codes', 'counter_values']:\n",
    "                            if hasattr(record_option.next_section.section, attr_name):\n",
    "                                attr_value = getattr(record_option.next_section.section, attr_name)\n",
    "                                detector_trig_data = detector_group_i.create_dataset(attr_name, data=attr_value)\n",
    "                                detector.append(detector_trig_data)\n",
    "                                if use_test_parse:\n",
    "                                    print(f'{attr_name}: {attr_value}')\n",
    "                        # Store each detector_trigger_rate in a higher level array\n",
    "                        detector_trig_array.append(detector)\n",
    "                        \n",
    "                    if type == 'veto_trigger_rates':\n",
    "                        # Store veto data in an array\n",
    "                        veto = []\n",
    "                        veto_group_i = veto_trig_grp.create_group(f'{type}_group_{veto_count}')\n",
    "                        veto_count += 1\n",
    "                        if use_test_parse:\n",
    "                            print(f'\\nVeto trigger rates data:')\n",
    "                        for attr_name in ['veto_trigger_header', 'len_to_next_header', 'clocking_interval',\n",
    "                                          'num_entries', 'detector_code', 'counter_value_det_code']:\n",
    "                            if hasattr(record_option.next_section.section, attr_name):\n",
    "                                attr_value = getattr(record_option.next_section.section, attr_name)\n",
    "                                veto_trig_data = veto_group_i.create_dataset(attr_name, data=attr_value)\n",
    "                                veto.append(veto_trig_data)\n",
    "                                if use_test_parse:\n",
    "                                    print(f'{attr_name}: {attr_value}')\n",
    "                        # Store each veto array in higher level array\n",
    "                        veto_trig_array.append(veto)\n",
    "\n",
    "        if use_test_parse:\n",
    "            # Print data about the arrays and their values\n",
    "            print(f'\\nDetector type set:     {detector_type_set}')\n",
    "            print(f'Unique detector codes: {len(detector_code_set)}')\n",
    "            print(f'Unique detector nums:  {len(detector_number_set)}')\n",
    "            print(f'Trace length set:      {trace_set}')\n",
    "            print(f'Charge length set:     {charge_set}')\n",
    "            print(f'Phonon length set:     {phonon_set}')\n",
    "            print(f'Error length set:      {error_set}')\n",
    "            print(f'Veto length set:       {veto_set}')\n",
    "            print(f'Num Events:            {len(event_set)}')\n",
    "            print(f'Charge array len:      {len(charge_pulse_array)}')\n",
    "            print(f'Phonon array len:      {len(phonon_pulse_array)}')\n",
    "            print(f'Veto array len:        {len(veto_array)}')\n",
    "            print(f'Error array len:       {len(error_array)}')\n",
    "\n",
    "\n",
    "#input_path  = \"../01120210_0727_F0114\"\n",
    "#input_path = \"/data3/afisher/test/01120210_0727_F0001\"\n",
    "#input_path = \"/data3/afisher/test/01130208_1838_F0006\"\n",
    "\n",
    "# This one should contain cut content\n",
    "input_path = 'minimal_sample.soudan'\n",
    "#output_path = \"/home/afisher@novateur.com/dataReaderWriter/NovateurData/01150212_1819_F0001_test_parse.hdf5\"\n",
    "output_path = \"minimal_sample.hdf5\"\n",
    "\n",
    "# For final files, save onto novateur network:\n",
    "#output_path = \"/data3/afisher/test/parsed_files/01150212_1819_F0001_parsed.hdf5\"\n",
    "\n",
    "parse_file(input_path, output_path, use_test_parse=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import h5py"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f = h5py.File(\"minimal_sample.hdf5\", 'r')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f.visititems(print_structure)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "myenv",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}