{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"The code in this notebook is based on the code as was originally written in [LSTM Preprocessing Point Picker](https://www.kaggle.com/code/seungmoklee/lstm-preprocessing-point-picker). The author of this notebook did a great job of setting a clear baseline.\n\nI modified the code in the following part:\n* Maximum pulse count is set to 96.\n* Remove the features r_err and z_err.\n* Remove all non-essential code and graphics. \n\nWith these few changes the output files only contain the features for the events as I use them in my [Tensorflow LSTM Model Training TPU](https://www.kaggle.com/code/rsmits/tensorflow-lstm-model-training-tpu) notebook and [Tensorflow LSTM Model Inference](https://www.kaggle.com/code/rsmits/tensorflow-lstm-model-inference) notebook.","metadata":{}},{"cell_type":"code","source":"# Data I/O and preprocessing\nimport numpy as np\nimport pandas as pd\nimport pyarrow.parquet as pq\n\n# System\nimport time\nimport os\nimport gc\nfrom tqdm.notebook import tqdm\n\n# multiprocessing\nimport multiprocessing","metadata":{"execution":{"iopub.status.busy":"2023-03-15T22:10:36.091173Z","iopub.execute_input":"2023-03-15T22:10:36.091750Z","iopub.status.idle":"2023-03-15T22:10:36.101913Z","shell.execute_reply.started":"2023-03-15T22:10:36.091715Z","shell.execute_reply":"2023-03-15T22:10:36.100463Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Data setting\ntrain_batch_id_first = 100\ntrain_batch_id_last = 101\ntrain_batch_ids = range(train_batch_id_first, train_batch_id_last + 1)\n\n# Feature Settings\nmax_pulse_count = 96\nn_features = 7  # time, charge, aux, x, y, z, rank \n\n# Directories\nhome_dir = \"/kaggle/input/icecube-neutrinos-in-deep-ice/\"\ntrain_format = home_dir + 'train/batch_{batch_id:d}.parquet'\npoint_picker_format = 'pp_mpc96_n7_batch_{batch_id:d}.npz'","metadata":{"execution":{"iopub.status.busy":"2023-03-15T22:10:36.105388Z","iopub.execute_input":"2023-03-15T22:10:36.105808Z","iopub.status.idle":"2023-03-15T22:10:36.127331Z","shell.execute_reply.started":"2023-03-15T22:10:36.105779Z","shell.execute_reply":"2023-03-15T22:10:36.125591Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Sensor Geometry Data\nsensor_geometry_df = pd.read_csv(home_dir + \"sensor_geometry.csv\")\n\n# X, Y, Z coordinates\nsensor_x = sensor_geometry_df.x\nsensor_y = sensor_geometry_df.y\nsensor_z = sensor_geometry_df.z\n\n# Detector constants\nc_const = 0.299792458  # speed of light [m/ns]\n\n# Min / Max information\nx_min = sensor_x.min()\nx_max = sensor_x.max()\ny_min = sensor_y.min()\ny_max = sensor_y.max()\nz_min = sensor_z.min()\nz_max = sensor_z.max()\n\n# Detector Valid Length\ndetector_length = np.sqrt((x_max - x_min)**2 + (y_max - y_min)**2 + (z_max - z_min)**2)\nt_valid_length = detector_length / c_const\n\nprint(f\"time valid length: {t_valid_length} ns\")","metadata":{"execution":{"iopub.status.busy":"2023-03-15T22:10:36.128840Z","iopub.execute_input":"2023-03-15T22:10:36.129875Z","iopub.status.idle":"2023-03-15T22:10:36.160276Z","shell.execute_reply.started":"2023-03-15T22:10:36.129838Z","shell.execute_reply":"2023-03-15T22:10:36.159298Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\n\n## Single event reader function\n\n- Pick-up important data points first\n    - Rank 3 (First)\n        - not aux, in valid time window\n    - Rank 2\n        - not aux, out of valid time window\n    - Rank 1\n        - aux, in valid time window\n    - Rank 0 (Last)\n        - aux, out of valid time window\n    - In each ranks, take pulses from highest charge\n\n\"\"\"\n\n# read single event from batch_meta_df\ndef read_event(event_idx, batch_meta_df, max_pulse_count, batch_df, train=True):\n    # read metadata\n    batch_id, first_pulse_index, last_pulse_index = batch_meta_df.iloc[event_idx][[\"batch_id\", \"first_pulse_index\", \"last_pulse_index\"]].astype(\"int\")\n\n    # read event\n    event_feature = batch_df[first_pulse_index:last_pulse_index + 1]\n    sensor_id = event_feature.sensor_id\n    \n    # merge features into single structured array\n    dtype = [(\"time\", \"float16\"),\n             (\"charge\", \"float16\"),\n             (\"auxiliary\", \"float16\"),\n             (\"x\", \"float16\"),\n             (\"y\", \"float16\"),\n             (\"z\", \"float16\"),\n             (\"rank\", \"short\")]\n    event_x = np.zeros(last_pulse_index - first_pulse_index + 1, dtype)\n\n    event_x[\"time\"] = event_feature.time.values - event_feature.time.min()\n    event_x[\"charge\"] = event_feature.charge.values\n    event_x[\"auxiliary\"] = event_feature.auxiliary.values\n\n    event_x[\"x\"] = sensor_geometry_df.x[sensor_id].values\n    event_x[\"y\"] = sensor_geometry_df.y[sensor_id].values\n    event_x[\"z\"] = sensor_geometry_df.z[sensor_id].values\n\n    # For long event, pick-up\n    if len(event_x) > max_pulse_count:\n        # Find valid time window\n        t_peak = event_x[\"time\"][event_x[\"charge\"].argmax()]\n        t_valid_min = t_peak - t_valid_length\n        t_valid_max = t_peak + t_valid_length\n\n        t_valid = (event_x[\"time\"] > t_valid_min) * (event_x[\"time\"] < t_valid_max)\n\n        # rank\n        event_x[\"rank\"] = 2 * (1 - event_x[\"auxiliary\"]) + (t_valid)\n\n        # sort by rank and charge (important goes to backward)\n        event_x = np.sort(event_x, order=[\"rank\", \"charge\"])\n\n        # pick-up from backward\n        event_x = event_x[-max_pulse_count:]\n\n        # resort by time\n        event_x = np.sort(event_x, order=\"time\")\n\n    # resort by time\n    event_x = np.sort(event_x, order=\"time\")\n        \n    # for train data, give angles together\n    azimuth, zenith = batch_meta_df.iloc[event_idx][[\"azimuth\", \"zenith\"]].astype(\"float16\")\n    event_y = np.array([azimuth, zenith], dtype=\"float16\")\n        \n    return event_idx, len(event_x), event_x, event_y","metadata":{"execution":{"iopub.status.busy":"2023-03-15T22:10:36.163001Z","iopub.execute_input":"2023-03-15T22:10:36.164933Z","iopub.status.idle":"2023-03-15T22:10:36.181216Z","shell.execute_reply.started":"2023-03-15T22:10:36.164878Z","shell.execute_reply":"2023-03-15T22:10:36.179171Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Read Train Meta Data\ntrain_meta_df = pd.read_parquet(home_dir + 'train_meta.parquet')\n\nbatch_counts = train_meta_df.batch_id.value_counts().sort_index()\n\nbatch_max_index = batch_counts.cumsum()\nbatch_max_index[train_meta_df.batch_id.min() - 1] = 0\nbatch_max_index = batch_max_index.sort_index()\n\ndef train_meta_df_spliter(batch_id):\n    return train_meta_df.loc[batch_max_index[batch_id - 1]:batch_max_index[batch_id] - 1]\n\nfor batch_id in train_batch_ids:\n    print(\"Reading batch \", batch_id, end=\"\")\n    # get batch meta data and data\n    batch_meta_df = train_meta_df_spliter(batch_id)\n    batch_df = pd.read_parquet(train_format.format(batch_id=batch_id))\n\n    # register pulses\n    batch_x = np.zeros((len(batch_meta_df), max_pulse_count, n_features), dtype=\"float16\")\n    batch_y = np.zeros((len(batch_meta_df), 2), dtype=\"float16\")\n    \n    batch_x[:, :, 2] = -1\n\n    def read_event_local(event_idx):\n        return read_event(event_idx, batch_meta_df, max_pulse_count, batch_df, train=True)\n\n    # Proces Events\n    iterator = range(len(batch_meta_df))\n    with multiprocessing.Pool() as pool:\n        for event_idx, pulse_count, event_x, event_y in pool.map(read_event_local, iterator):\n            batch_x[event_idx, :pulse_count, 0] = event_x[\"time\"]\n            batch_x[event_idx, :pulse_count, 1] = event_x[\"charge\"]\n            batch_x[event_idx, :pulse_count, 2] = event_x[\"auxiliary\"]\n            batch_x[event_idx, :pulse_count, 3] = event_x[\"x\"]\n            batch_x[event_idx, :pulse_count, 4] = event_x[\"y\"]\n            batch_x[event_idx, :pulse_count, 5] = event_x[\"z\"]\n            batch_x[event_idx, :pulse_count, 6] = event_x[\"rank\"]\n\n            batch_y[event_idx] = event_y\n\n    del batch_meta_df, batch_df\n    \n    # Save\n    print(\" DONE! Saving...\")\n    np.savez(point_picker_format.format(batch_id=batch_id), x=batch_x, y=batch_y)","metadata":{"execution":{"iopub.status.busy":"2023-03-15T22:10:36.183327Z","iopub.execute_input":"2023-03-15T22:10:36.183833Z"},"trusted":true},"execution_count":null,"outputs":[]}]}