{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.12.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceType":"competition","sourceId":101849,"databundleVersionId":13093295}],"dockerImageVersionId":31328,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\n\n# 1. Load and reshape the FGS1 signal\n# Spatial dimensions: 32x32\nfgs1_path = \"/kaggle/input/competitions/ariel-data-challenge-2025/train/1010375142/FGS1_signal_0.parquet\"\nfgs1_df = pd.read_parquet(fgs1_path)\nfgs1_images = fgs1_df.to_numpy().reshape(-1, 32, 32)\n\nprint(f\"FGS1 3D Cube Shape: {fgs1_images.shape}\") \n# Expected output: (135000, 32, 32)\n\n# 2. Load and reshape the AIRS-CH0 signal\n# Spatial dimensions: 32x356\nairs_path = \"/kaggle/input/competitions/ariel-data-challenge-2025/train/1010375142/AIRS-CH0_signal_0.parquet\"\nairs_df = pd.read_parquet(airs_path)\nairs_images = airs_df.to_numpy().reshape(-1, 32, 356)\n\nprint(f\"AIRS-CH0 3D Cube Shape: {airs_images.shape}\") \n# Expected output: (11250, 32, 356)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:09.521766Z","iopub.execute_input":"2026-04-25T16:37:09.522089Z","iopub.status.idle":"2026-04-25T16:37:15.202681Z","shell.execute_reply.started":"2026-04-25T16:37:09.522056Z","shell.execute_reply":"2026-04-25T16:37:15.201589Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# **Exploratory Data Analysis (EDA)**","metadata":{}},{"cell_type":"code","source":"import matplotlib.pyplot as plt\n\n# --- EDA 1: Visualize a Single Sensor Frame ---\n# Let's look at the very first snapshot captured by the AIRS-CH0 instrument\nplt.figure(figsize=(12, 4))\nplt.imshow(fgs1_images[0], cmap='viridis', aspect='auto')\nplt.colorbar(label='Signal Intensity')\nplt.title('fgs1 2D Sensor Frame (Time Step 0)')\nplt.xlabel('Spectral Dimension (Wavelengths)')\nplt.ylabel('Spatial Dimension')\nplt.show()\n\n# --- EDA 1: Visualize a Single Sensor Frame ---\n# Let's look at the very first snapshot captured by the AIRS-CH0 instrument\nplt.figure(figsize=(12, 4))\nplt.imshow(airs_images[0], cmap='viridis', aspect='auto')\nplt.colorbar(label='Signal Intensity')\nplt.title('AIRS-CH0 2D Sensor Frame (Time Step 0)')\nplt.xlabel('Spectral Dimension (Wavelengths)')\nplt.ylabel('Spatial Dimension')\nplt.show()\n\n# --- EDA 2: Plot the Light Curve (Transit Dip) ---\n# To see the planet pass in front of the star, we average the pixels across the spatial dimensions.\n# This gives us the total brightness for each time step.\nfgs1_lightcurve = fgs1_images.mean(axis=(1, 2))\n\n# Note: The raw signal is extremely noisy. In practice, you will want to apply a rolling mean.\nfgs1_smoothed = pd.Series(fgs1_lightcurve).rolling(window=500, center=True).mean()\n\nplt.figure(figsize=(12, 5))\nplt.plot(fgs1_lightcurve, alpha=0.3, color='gray', label='Raw Signal')\nplt.plot(fgs1_smoothed, color='blue', linewidth=2, label='Smoothed Light Curve')\nplt.title('FGS1 Light Curve: Spotting the Exoplanet Transit')\nplt.xlabel('Time Step')\nplt.ylabel('Mean Pixel Intensity')\nplt.legend()\nplt.show()\n\n# --- EDA 2: Plot the Light Curve (Transit Dip) ---\n# To see the planet pass in front of the star, we average the pixels across the spatial dimensions.\n# This gives us the total brightness for each time step.\nairs_lightcurve = airs_images.mean(axis=(1, 2))\n\n# Note: The raw signal is extremely noisy. In practice, you will want to apply a rolling mean.\nairs_smoothed = pd.Series(airs_lightcurve).rolling(window=500, center=True).mean()\n\nplt.figure(figsize=(12, 5))\nplt.plot(airs_lightcurve, alpha=0.3, color='gray', label='Raw Signal')\nplt.plot(airs_smoothed, color='blue', linewidth=2, label='Smoothed Light Curve')\nplt.title('AIRS Light Curve: Spotting the Exoplanet Transit')\nplt.xlabel('Time Step')\nplt.ylabel('Mean Pixel Intensity')\nplt.legend()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:15.20413Z","iopub.execute_input":"2026-04-25T16:37:15.20456Z","iopub.status.idle":"2026-04-25T16:37:20.282923Z","shell.execute_reply.started":"2026-04-25T16:37:15.204514Z","shell.execute_reply":"2026-04-25T16:37:20.281775Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# **Caliberation**","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport itertools\nfrom astropy.stats import sigma_clip\n\n# --- 1. Analog-to-Digital Conversion ---\ndef ADC_convert(signal, gain=0.4369, offset=-1000):\n    signal = signal.astype(np.float64)\n    signal /= gain\n    signal += offset\n    return signal\n\n# --- 2. Mask Hot/Dead Pixels ---\ndef mask_hot_dead(signal, dead, dark):\n    # Sigma clipping to find hot pixels in the dark frame\n    hot = sigma_clip(dark, sigma=5, maxiters=5).mask\n    \n    # Broadcast masks to match the time-series shape\n    hot = np.tile(hot, (signal.shape[0], 1, 1))\n    dead = np.tile(dead, (signal.shape[0], 1, 1))\n    \n    # Apply masks\n    signal = np.ma.masked_where(dead, signal)\n    signal = np.ma.masked_where(hot, signal)\n    return signal\n\n# --- 3. Linearity Correction ---\ndef apply_linear_corr(linear_corr, clean_signal):\n    linear_corr = np.flip(linear_corr, axis=0)\n    for x, y in itertools.product(range(clean_signal.shape[1]), range(clean_signal.shape[2])):\n        poli = np.poly1d(linear_corr[:, x, y])\n        clean_signal[:, x, y] = poli(clean_signal[:, x, y])\n    return clean_signal\n\n# --- 4. Dark Current Subtraction ---\ndef clean_dark(signal, dead, dark, dt):\n    dark = np.ma.masked_where(dead, dark)\n    dark = np.tile(dark, (signal.shape[0], 1, 1))\n    # Multiply dark current by integration time (dt)\n    signal -= dark * dt[:, np.newaxis, np.newaxis]\n    return signal\n\n# --- 5. Correlated Double Sampling (CDS) ---\ndef get_cds(signal):\n    # End of exposure (odds) minus Start of exposure (evens)\n    cds = signal[1::2, :, :] - signal[::2, :, :]\n    return cds\n\n# --- 6. Flat Field Correction ---\ndef correct_flat_field(flat, dead, signal):\n    flat = np.ma.masked_where(dead, flat)\n    flat = np.tile(flat, (signal.shape[0], 1, 1))\n    signal = signal / flat\n    return signal\n\n# ==========================================\n# MASTER CALIBRATION FUNCTION\n# ==========================================\ndef calibrate_signal(raw_signal, dt_array, cal_dict, config):\n    \"\"\"\n    Applies the full calibration sequence to a single raw signal array.\n    \n    Parameters:\n    - raw_signal: 3D numpy array (time_steps, height, width)\n    - dt_array: 1D numpy array of integration times matching time_steps\n    - cal_dict: Dictionary containing 'dark', 'dead', 'flat', 'linear_corr' arrays\n    - config: Dictionary of boolean flags to toggle steps (e.g., {'DO_MASK': True, ...})\n    \"\"\"\n    \n    # Step 1: ADC Conversion\n    signal = ADC_convert(raw_signal)\n    \n    # Step 2: Masking\n    if config.get('DO_MASK', True):\n        signal = mask_hot_dead(signal, cal_dict['dead'], cal_dict['dark'])\n        \n    # Step 3: Linearity Correction\n    if config.get('DO_THE_NL_CORR', False):\n        signal = apply_linear_corr(cal_dict['linear_corr'], signal)\n        \n    # Step 4: Dark Current\n    if config.get('DO_DARK', True):\n        signal = clean_dark(signal, cal_dict['dead'], cal_dict['dark'], dt_array)\n        \n    # Step 5: CDS (Transforms shape from N to N/2)\n    signal_cds = get_cds(signal)\n    \n    # Step 6: Flat Field\n    if config.get('DO_FLAT', True):\n        signal_cds = correct_flat_field(cal_dict['flat'], cal_dict['dead'], signal_cds)\n        \n    return signal_cds","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:20.284376Z","iopub.execute_input":"2026-04-25T16:37:20.284853Z","iopub.status.idle":"2026-04-25T16:37:20.983999Z","shell.execute_reply.started":"2026-04-25T16:37:20.284819Z","shell.execute_reply":"2026-04-25T16:37:20.983001Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\n\n# ==========================================\n# 1. LOAD DATA INTO MEMORY\n# ==========================================\n# Define paths\nobs_id = \"1010375142\"\n\nfgs1_signal_path = f\"/kaggle/input/competitions/ariel-data-challenge-2025/train/{obs_id}/FGS1_signal_0.parquet\"\nfgs1_cal_dir = f\"/kaggle/input/competitions/ariel-data-challenge-2025/train/{obs_id}/FGS1_calibration_0\"\n\n# Load raw signal (reshape to 3D: time, height, width)\nfgs1_raw = pd.read_parquet(fgs1_signal_path).values.astype(np.float64).reshape(-1, 32, 32)\n\n# Load calibration frames (reshape to 2D: height, width)\ndark = pd.read_parquet(fgs1_cal_dir + \"/dark.parquet\").values.astype(np.float64).reshape(32, 32)\ndead = pd.read_parquet(fgs1_cal_dir + \"/dead.parquet\").values.astype(np.float64).reshape(32, 32)\nflat = pd.read_parquet(fgs1_cal_dir + \"/flat.parquet\").values.astype(np.float64).reshape(32, 32)\n\n# Recreate the integration time array (dt) based on your original script\ndt_fgs1 = np.ones(len(fgs1_raw)) * 0.1\ndt_fgs1[1::2] += 0.1 \n\n# ==========================================\n# 2. APPLY CALIBRATION PIPELINE\n# ==========================================\nprint(f\"Original shape: {fgs1_raw.shape}\")\n\n# Step 1: Convert Analog to Digital\nclean_fgs1 = ADC_convert(fgs1_raw)\n\n# Step 2: Mask Dead & Hot Pixels\nclean_fgs1 = mask_hot_dead(clean_fgs1, dead, dark)\n\n# Step 3: Subtract Dark Current (Thermal Noise)\nclean_fgs1 = clean_dark(clean_fgs1, dead, dark, dt_fgs1)\n\n# Step 4: Correlated Double Sampling (CDS)\n# NOTE: This subtracts pairs of frames, so it will cut your time steps exactly in half!\nclean_fgs1 = get_cds(clean_fgs1)\n\n# Step 5: Flat Field Correction (Pixel sensitivity)\nclean_fgs1 = correct_flat_field(flat, dead, clean_fgs1)\n\nprint(f\"Calibrated shape (after CDS): {clean_fgs1.shape}\")\n\n\n# ==========================================\n# 3. EDA PLOTS (CALIBRATED DATA)\n# ==========================================\n\n# --- EDA 1: Visualize the Calibrated 2D Sensor Frame ---\nplt.figure(figsize=(6, 5))\n# We use [0] to look at the very first cleaned frame\nplt.imshow(clean_fgs1[0], cmap='viridis', aspect='auto')\nplt.colorbar(label='Calibrated Signal Intensity')\nplt.title('FGS1 Calibrated 2D Sensor Frame (Time Step 0)')\nplt.xlabel('Pixel Width')\nplt.ylabel('Pixel Height')\nplt.show()\n\n# --- EDA 2: Plot the Calibrated Light Curve ---\n# Average the pixels across the spatial dimensions (axis 1 and 2)\nfgs1_calibrated_lightcurve = clean_fgs1.mean(axis=(1, 2))\n\n# Apply a rolling mean to smooth out remaining noise. \n# Window is 250 instead of 500 because CDS cut the frame count in half.\nfgs1_smoothed = pd.Series(fgs1_calibrated_lightcurve).rolling(window=250, center=True).mean()\n\nplt.figure(figsize=(12, 5))\nplt.plot(fgs1_calibrated_lightcurve, alpha=0.3, color='gray', label='Calibrated Signal (CDS)')\nplt.plot(fgs1_smoothed, color='red', linewidth=2, label='Smoothed Calibrated Light Curve')\nplt.title('FGS1 CALIBRATED Light Curve: Exoplanet Transit')\nplt.xlabel('Time Step (CDS Frames)')\nplt.ylabel('Mean Pixel Intensity')\nplt.legend()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:20.986299Z","iopub.execute_input":"2026-04-25T16:37:20.986834Z","iopub.status.idle":"2026-04-25T16:37:31.892526Z","shell.execute_reply.started":"2026-04-25T16:37:20.986795Z","shell.execute_reply":"2026-04-25T16:37:31.891388Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport os\n\n# ==========================================\n# 1. LOAD DATA & CHOP EDGES\n# ==========================================\nobs_id = \"1010375142\"\npath_folder = \"/kaggle/input/competitions/ariel-data-challenge-2025/\" # Update if needed\n\nairs_signal_path = f\"{path_folder}train/{obs_id}/AIRS-CH0_signal_0.parquet\"\nairs_cal_dir = f\"{path_folder}train/{obs_id}/AIRS-CH0_calibration_0\"\n\n# Set the limits to chop the wavelength dimension (from your original code)\ncut_inf, cut_sup = 39, 321\n\n# Load and reshape raw signal, then immediately slice the wavelength axis\nairs_raw = pd.read_parquet(airs_signal_path).values.astype(np.float64).reshape(-1, 32, 356)\nairs_raw = airs_raw[:, :, cut_inf:cut_sup]\n\n# Load and slice calibration frames to match the new 32x282 shape\ndark = pd.read_parquet(airs_cal_dir + \"/dark.parquet\").values.astype(np.float64).reshape(32, 356)[:, cut_inf:cut_sup]\ndead = pd.read_parquet(airs_cal_dir + \"/dead.parquet\").values.astype(np.float64).reshape(32, 356)[:, cut_inf:cut_sup]\nflat = pd.read_parquet(airs_cal_dir + \"/flat.parquet\").values.astype(np.float64).reshape(32, 356)[:, cut_inf:cut_sup]\n\n# Load integration time from axis_info\naxis_info = pd.read_parquet(os.path.join(path_folder, 'axis_info.parquet'))\ndt_airs = axis_info['AIRS-CH0-integration_time'].dropna().values\n# Note: For a single observation chunk, dt_airs might need to be trimmed to match the frame count (len(airs_raw))\n# For this example, we'll assume it matches or create a dummy if axis_info isn't loaded properly:\nif len(dt_airs) < len(airs_raw):\n    # Fallback just in case axis_info doesn't map 1:1 in your current notebook state\n    dt_airs = np.ones(len(airs_raw)) * 0.1 \ndt_airs[1::2] += 0.1\n\n# ==========================================\n# 2. APPLY CALIBRATION PIPELINE\n# ==========================================\nprint(f\"Original chopped shape: {airs_raw.shape}\")\n\n# Run through the pipeline steps\nclean_airs = ADC_convert(airs_raw)\nclean_airs = mask_hot_dead(clean_airs, dead, dark)\nclean_airs = clean_dark(clean_airs, dead, dark, dt_airs)\nclean_airs = get_cds(clean_airs)\nclean_airs = correct_flat_field(flat, dead, clean_airs)\n\nprint(f\"Calibrated shape (after CDS): {clean_airs.shape}\")\n\n# ==========================================\n# 3. EDA PLOTS (CALIBRATED DATA)\n# ==========================================\n\n# --- EDA 1: Visualize the Calibrated 2D Sensor Frame ---\nplt.figure(figsize=(10, 4))\n# We use vmin and vmax to stretch the contrast, making the spectrum easier to see\nplt.imshow(clean_airs[0], cmap='viridis', aspect='auto', vmin=np.percentile(clean_airs[0], 5), vmax=np.percentile(clean_airs[0], 95))\nplt.colorbar(label='Calibrated Signal Intensity')\nplt.title('AIRS-CH0 Calibrated Spectrum Frame (Time Step 0)')\nplt.xlabel('Wavelength Pixel (Trimmed)')\nplt.ylabel('Spatial Pixel')\nplt.show()\n\n# --- EDA 2: Plot the Calibrated Light Curve (Broadband) ---\n# To get the \"white light\" (broadband) curve, we sum/average across ALL wavelengths\nairs_calibrated_lightcurve = clean_airs.mean(axis=(1, 2))\n\n# Smooth the light curve (AIRS has fewer frames than FGS1, so use a smaller window)\nairs_smoothed = pd.Series(airs_calibrated_lightcurve).rolling(window=50, center=True).mean()\n\nplt.figure(figsize=(12, 5))\nplt.plot(airs_calibrated_lightcurve, alpha=0.4, color='gray', label='Calibrated Signal (CDS)')\nplt.plot(airs_smoothed, color='purple', linewidth=2, label='Smoothed Broadband Light Curve')\nplt.title('AIRS-CH0 Calibrated Broadband Light Curve')\nplt.xlabel('Time Step (CDS Frames)')\nplt.ylabel('Mean Pixel Intensity')\nplt.legend()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:31.894229Z","iopub.execute_input":"2026-04-25T16:37:31.894577Z","iopub.status.idle":"2026-04-25T16:37:39.752271Z","shell.execute_reply.started":"2026-04-25T16:37:31.894544Z","shell.execute_reply":"2026-04-25T16:37:39.751198Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# **Plot the OUTPUT Transmission Spectrum**","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport matplotlib.pyplot as plt\n\n# 1. Load the labels from train.csv (update path to match your Kaggle directory)\nlabels_path = \"/kaggle/input/competitions/ariel-data-challenge-2025/train.csv\"\ntrain_df = pd.read_csv(labels_path)\n\n# 2. Extract the targets for observation 1010375142\nobs_id = 1010375142\n\n# Filter the dataframe for our specific planet ID\nplanet_data = train_df[train_df['planet_id'] == obs_id]\n\n# Drop the planet_id column so we are left with ONLY the 283 target values\nplanet_targets = planet_data.drop(columns=['planet_id']).values[0]\n\n# Ensure we have exactly 283 values\nprint(f\"Number of targets extracted: {len(planet_targets)}\")\n\n# 3. Separate AIRS-CH0 (Infrared) and FGS1 (Visible/Broadband)\nairs_targets = planet_targets[:282]\nfgs1_target = planet_targets[282]\n\n# 4. Plot the Spectrum\nplt.figure(figsize=(14, 6))\n\nplt.plot(range(1, 283), airs_targets, marker='o', markersize=4, linestyle='-', color='purple', alpha=0.7, label='AIRS-CH0 (Infrared Wavelengths)')\nplt.scatter([283], [fgs1_target], color='red', marker='*', s=200, label='FGS1 (Visible Anchor Point)')\n\nplt.title(f'Exoplanet Transmission Spectrum (Obs ID: {obs_id})', fontsize=16)\nplt.xlabel('Wavelength Target Bins (1-282 = Infrared, 283 = Visible)', fontsize=12)\nplt.ylabel('Transit Depth (Relative Planet Size)', fontsize=12)\nplt.grid(True, alpha=0.3)\nplt.legend(fontsize=12)\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:39.75357Z","iopub.execute_input":"2026-04-25T16:37:39.753984Z","iopub.status.idle":"2026-04-25T16:37:40.137577Z","shell.execute_reply.started":"2026-04-25T16:37:39.753954Z","shell.execute_reply":"2026-04-25T16:37:40.136696Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport matplotlib.pyplot as plt\nimport numpy as np\n\n# 1. Load the labels from train.csv\nlabels_path = \"/kaggle/input/competitions/ariel-data-challenge-2025/train.csv\"\ntrain_df = pd.read_csv(labels_path)\n\n# 2. Extract the targets for observation 1010375142\nobs_id = 1010375142\nplanet_data = train_df[train_df['planet_id'] == obs_id]\nplanet_targets = planet_data.drop(columns=['planet_id']).values[0]\n\nairs_targets = planet_targets[:282]\nfgs1_target = planet_targets[282]\n\n# --- NEW: CALCULATE THE BASELINE ---\n# We use the 5th percentile of the infrared spectrum to represent the \"bare planet\" or \"cloud deck\" floor\nbaseline_depth = np.percentile(airs_targets, 5)\n\n# 4. Plot the Spectrum\nplt.figure(figsize=(14, 6))\n\n# Plot the spectrum and FGS1 point\nplt.plot(range(1, 283), airs_targets, marker='o', markersize=4, linestyle='-', color='purple', alpha=0.7, label='AIRS-CH0 (Atmosphere + Core)')\nplt.scatter([283], [fgs1_target], color='red', marker='*', s=200, label='FGS1 (Visible Anchor)')\n\n# --- NEW: DRAW THE BASELINE ---\nplt.axhline(y=baseline_depth, color='black', linestyle='--', linewidth=2, label=f'Baseline/Core Depth (~{baseline_depth:.4f})')\n\n# Fill the area between the baseline and the peaks to highlight the molecules!\nplt.fill_between(range(1, 283), baseline_depth, airs_targets, where=(airs_targets > baseline_depth), color='purple', alpha=0.1, label='Molecular Absorption')\n\nplt.title(f'Exoplanet Transmission Spectrum (Obs ID: {obs_id})', fontsize=16)\nplt.xlabel('Wavelength Target Bins (1-282 = Infrared, 283 = Visible)', fontsize=12)\nplt.ylabel('Transit Depth (Relative Planet Size)', fontsize=12)\nplt.grid(True, alpha=0.3)\nplt.legend(fontsize=12)\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:40.138533Z","iopub.execute_input":"2026-04-25T16:37:40.13878Z","iopub.status.idle":"2026-04-25T16:37:40.583479Z","shell.execute_reply.started":"2026-04-25T16:37:40.138755Z","shell.execute_reply":"2026-04-25T16:37:40.582634Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**EDA 1: Raw vs. Calibrated (Side-by-Side)**","metadata":{}},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport numpy as np\n\n# Choose a frame to look at (e.g., the 10th time step)\nframe_idx = 10\n\nfig, axes = plt.subplots(1, 2, figsize=(15, 5))\n\n# Plot 1: Raw, uncalibrated data\n# We use vmin/vmax to handle extreme outlier pixels so the image doesn't just look black\nim1 = axes[0].imshow(airs_raw[frame_idx], cmap='magma', aspect='auto', \n                     vmin=np.percentile(airs_raw[frame_idx], 5), \n                     vmax=np.percentile(airs_raw[frame_idx], 95))\naxes[0].set_title('RAW AIRS-CH0 Image (Frame 10)')\naxes[0].set_xlabel('Wavelength Pixels')\naxes[0].set_ylabel('Spatial Pixels')\nfig.colorbar(im1, ax=axes[0], label='Analog Signal')\n\n# Plot 2: Cleaned, calibrated data\nim2 = axes[1].imshow(clean_airs[frame_idx], cmap='magma', aspect='auto',\n                     vmin=np.percentile(clean_airs[frame_idx], 5), \n                     vmax=np.percentile(clean_airs[frame_idx], 95))\naxes[1].set_title('CALIBRATED AIRS-CH0 Image (Frame 10)')\naxes[1].set_xlabel('Wavelength Pixels')\naxes[1].set_ylabel('Spatial Pixels')\nfig.colorbar(im2, ax=axes[1], label='Calibrated Intensity')\n\nplt.tight_layout()\nplt.show()\n\n\n\n\nimport matplotlib.pyplot as plt\nimport numpy as np\n\n# Choose a frame to look at (e.g., the 10th time step)\nframe_idx = 10\n\nfig, axes = plt.subplots(1, 2, figsize=(15, 5))\n\n# Plot 1: Raw, uncalibrated data\n# We use vmin/vmax to handle extreme outlier pixels so the image doesn't just look black\nim1 = axes[0].imshow(fgs1_raw[frame_idx], cmap='magma', aspect='auto', \n                     vmin=np.percentile(fgs1_raw[frame_idx], 5), \n                     vmax=np.percentile(fgs1_raw[frame_idx], 95))\naxes[0].set_title('RAW FGS1 Image (Frame 10)')\naxes[0].set_xlabel('Wavelength Pixels')\naxes[0].set_ylabel('Spatial Pixels')\nfig.colorbar(im1, ax=axes[0], label='Analog Signal')\n\n# Plot 2: Cleaned, calibrated data\nim2 = axes[1].imshow(clean_fgs1[frame_idx], cmap='magma', aspect='auto',\n                     vmin=np.percentile(clean_fgs1[frame_idx], 5), \n                     vmax=np.percentile(clean_fgs1[frame_idx], 95))\naxes[1].set_title('CALIBRATED FGS1 Image (Frame 10)')\naxes[1].set_xlabel('Wavelength Pixels')\naxes[1].set_ylabel('Spatial Pixels')\nfig.colorbar(im2, ax=axes[1], label='Calibrated Intensity')\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T17:00:16.406346Z","iopub.execute_input":"2026-04-25T17:00:16.406743Z","iopub.status.idle":"2026-04-25T17:00:17.569878Z","shell.execute_reply.started":"2026-04-25T17:00:16.406707Z","shell.execute_reply":"2026-04-25T17:00:17.568731Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**EDA 2: The \"Waterfall\" Plot (Time vs. Wavelength)**","metadata":{}},{"cell_type":"code","source":"# 1. Compress the spatial dimension (axis=1) by summing the pixels vertically\n# Shape goes from (Time, 32, Wavelengths) -> (Time, Wavelengths)\nwaterfall_data = clean_airs.sum(axis=1)\n\n# 2. Plot the Waterfall Spectrogram\nplt.figure(figsize=(12, 8))\n\n# We plot the 2D array. \nplt.imshow(waterfall_data, cmap='viridis', aspect='auto', origin='lower')\n\nplt.colorbar(label='Total Intensity (Summed across spatial pixels)')\nplt.title('AIRS-CH0 Spectrogram: Time vs. Wavelength', fontsize=16)\nplt.xlabel('Wavelength Pixel Bin', fontsize=12)\nplt.ylabel('Time Step (CDS Frames)', fontsize=12)\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:41.945788Z","iopub.execute_input":"2026-04-25T16:37:41.946131Z","iopub.status.idle":"2026-04-25T16:37:43.015757Z","shell.execute_reply.started":"2026-04-25T16:37:41.946099Z","shell.execute_reply":"2026-04-25T16:37:43.014975Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# **EDA** ","metadata":{}},{"cell_type":"markdown","source":"**1. Understanding the Dataset First**\n* Each timestamp represents a single capture event (a frame). The sensor opens its electronic shutter, collects photons (light) for a specific fraction of a second, saves that 2D array of pixels, and then resets.\n* Single camera/source or multiple?\nMultiple. As we discussed, there are two independent cameras (FGS1 and AIRS-CH0) taking pictures at the exact same time, but at completely different frame rates.\n* not temporal labels (like \"anomaly at frame 500\"). The label is a single, static 283-number array (the transmission spectrum) for the entire observation sequence, located in train.csv.","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport os\n\ndef run_temporal_diagnostics(dt_gaps, instrument_name):\n    \"\"\"\n    Runs a temporal consistency check on an array of time gaps (integration times).\n    \"\"\"\n    print(f\"========== TEMPORAL REPORT: {instrument_name} ==========\")\n    \n    # Calculate absolute timestamps by cumulatively summing the gaps\n    timestamps = np.cumsum(dt_gaps)\n    \n    # ---------------------------------------------------------\n    # 1. UNDERSTAND THE DATASET\n    # ---------------------------------------------------------\n    total_frames = len(dt_gaps)\n    unique_gaps = np.unique(np.round(dt_gaps, 4)) # Rounding to handle floating point errors\n    \n    print(\"1. Dataset Overview:\")\n    print(f\"   -> Total capture events (frames): {total_frames}\")\n    \n    if len(unique_gaps) == 1:\n        print(f\"   -> Spacing: EVENLY spaced. Every gap is exactly {unique_gaps[0]}s\")\n    else:\n        print(f\"   -> Spacing: IRREGULAR. Detected gap sizes (seconds): {unique_gaps}\")\n        \n    print(\"   -> Source: Single instrument view (Part of a multi-camera system)\")\n    print(\"   -> Labels Available?: No frame-by-frame labels. One global label array per observation.\\n\")\n\n    # ---------------------------------------------------------\n    # 2. CHECK TEMPORAL CONSISTENCY\n    # ---------------------------------------------------------\n    print(\"2. Temporal Consistency Check:\")\n    \n    # Check for missing gaps (NaNs or 0s)\n    missing_count = np.isnan(dt_gaps).sum() + (dt_gaps == 0).sum()\n    print(f\"   -> Missing/Zero gaps (Dropped frames): {missing_count}\")\n    \n    # Check for duplicate absolute timestamps \n    # (If time doesn't move forward, the sensor froze)\n    duplicate_timestamps = len(timestamps) - len(np.unique(timestamps))\n    print(f\"   -> Duplicate timestamps (Sensor freezes): {duplicate_timestamps}\")\n    \n    # Check for negative time leaps\n    negative_gaps = (dt_gaps < 0).sum()\n    print(f\"   -> Negative time leaps (Time travel bugs): {negative_gaps}\\n\")\n\n    # ---------------------------------------------------------\n    # 3. PLOT: FRAME INDEX VS TIMESTAMP GAPS\n    # ---------------------------------------------------------\n    plt.figure(figsize=(12, 4))\n    \n    # Plotting only a slice (e.g., first 100 frames) so the pattern isn't compressed into a blur\n    frames_to_plot = min(100, total_frames)\n    frame_indices = np.arange(frames_to_plot)\n    \n    plt.plot(frame_indices, dt_gaps[:frames_to_plot], marker='o', linestyle='-', color='teal')\n    \n    plt.title(f'Temporal Gaps vs Frame Index: {instrument_name} (First {frames_to_plot} Frames)', fontsize=14)\n    plt.xlabel('Frame Index', fontsize=12)\n    plt.ylabel('Gap to Next Frame (seconds)', fontsize=12)\n    \n    # Add a buffer to the y-axis limits to clearly see the upper and lower bounds\n    plt.ylim(min(dt_gaps) - 0.05, max(dt_gaps) + 0.05)\n    plt.grid(True, alpha=0.4)\n    plt.tight_layout()\n    plt.show()\n\n# ==========================================\n# EXECUTE THE SCRIPT ON KAGGLE DATA\n# ==========================================\npath_folder = \"/kaggle/input/competitions/ariel-data-challenge-2025/\"\n\n# Load the axis_info file to get the exact timing arrays\naxis_info_path = os.path.join(path_folder, 'axis_info.parquet')\n\nif os.path.exists(axis_info_path):\n    axis_info = pd.read_parquet(axis_info_path)\n    \n    # Extract the raw gaps for AIRS-CH0\n    # Note: We simulate the raw sensor readout overhead (adding 0.1s to every alternate frame)\n    airs_base_gaps = axis_info['AIRS-CH0-integration_time'].dropna().values\n    airs_raw_gaps = np.copy(airs_base_gaps)\n    airs_raw_gaps[1::2] += 0.1 \n    \n    # Run the diagnostics!\n    run_temporal_diagnostics(airs_raw_gaps, instrument_name=\"AIRS-CH0\")\nelse:\n    print(f\"File not found at {axis_info_path}. Please check your Kaggle input path.\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:43.018397Z","iopub.execute_input":"2026-04-25T16:37:43.018871Z","iopub.status.idle":"2026-04-25T16:37:43.21302Z","shell.execute_reply.started":"2026-04-25T16:37:43.018841Z","shell.execute_reply":"2026-04-25T16:37:43.212021Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**The Visual Inspection Code**","metadata":{}},{"cell_type":"markdown","source":"* The Star's Light: The thick, glowing yellow/purple horizontal band is the raw starlight smeared across different wavelengths.\n\n* Uncleaned Hardware Defects: The sharp white needles above and below the band are \"hot pixels\" (permanently stuck-on sensor pixels). Because they appear in every frame, it proves your current calibration pipeline did not fully clean the noise.\n\n* Good Stability: The glowing band does not bounce up and down between consecutive frames, meaning the telescope was pointing steadily without much vibration (jitter).\n\n* The Invisible Transit: You cannot see the light dimming from the first frame to the middle frame because the exoplanet only blocks about 0.1% of the light—impossible for human eyes to spot without math.","metadata":{}},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport numpy as np\nimport random\n\n# Assuming clean_airs is your 3D array of shape (Time, 32, 356)\n# For testing purposes, if you don't have it loaded, we will define the total frames\ntotal_frames = clean_airs.shape[0]\n\n# --- Setup indices for inspection ---\n\n# 1. Consecutive Frames (e.g., frames 500, 501, 502)\nstart_consec = 500\nconsec_indices = [start_consec, start_consec + 1, start_consec + 2]\n\n# 2. Random Samples\nrandom_indices = random.sample(range(100, total_frames - 100), 3)\nrandom_indices.sort()\n\n# 3. Edge Cases (First, Last, and potentially anomalous frame)\nfirst_frame = 0\nlast_frame = total_frames - 1\n# Let's find an \"anomalous\" frame by looking for the one with the highest maximum pixel value \n# (often indicative of a cosmic ray hit or a severe jitter spike)\nmax_pixels_per_frame = np.max(clean_airs, axis=(1, 2))\nanomaly_frame = np.argmax(max_pixels_per_frame)\n\nedge_indices = [first_frame, anomaly_frame, last_frame]\n\n# --- Helper Function to Plot a Row of Images ---\ndef plot_image_row(indices, title, cmap='magma'):\n    fig, axes = plt.subplots(1, 3, figsize=(15, 4))\n    fig.suptitle(title, fontsize=14, fontweight='bold')\n    \n    for i, frame_idx in enumerate(indices):\n        # Calculate percentiles for vmin/vmax to handle extreme outliers \n        vmin = np.percentile(clean_airs[frame_idx], 2)\n        vmax = np.percentile(clean_airs[frame_idx], 98)\n        \n        im = axes[i].imshow(clean_airs[frame_idx], cmap=cmap, aspect='auto', vmin=vmin, vmax=vmax)\n        axes[i].set_title(f\"Frame {frame_idx}\")\n        axes[i].set_xlabel(\"Wavelength Pixels\")\n        if i == 0:\n            axes[i].set_ylabel(\"Spatial Pixels\")\n            \n    # Add a colorbar to the last plot to see intensity scales\n    fig.colorbar(im, ax=axes.ravel().tolist(), label=\"Pixel Intensity\")\n    plt.show()\n\n# --- Execute the Plots ---\nprint(\"Executing Visual Inspection...\")\n\nplot_image_row(consec_indices, \"1. Consecutive Frames (Look for smooth continuity or sudden shifts)\")\nplot_image_row(random_indices, \"2. Random Samples Across Time (Look for gradual lighting changes)\")\nplot_image_row(edge_indices, \"3. Edge Cases: First, Max Intensity (Anomaly), Last\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:43.214245Z","iopub.execute_input":"2026-04-25T16:37:43.214539Z","iopub.status.idle":"2026-04-25T16:37:44.889419Z","shell.execute_reply.started":"2026-04-25T16:37:43.214506Z","shell.execute_reply":"2026-04-25T16:37:44.888283Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":" **The Image Statistics Code**","metadata":{}},{"cell_type":"markdown","source":"* The Exoplanet's Shadow (Blue Graph): The large, smooth \"U-shape\" dip in the Mean Intensity plot is the actual transit. You are seeing the average light drop as the planet moves in front of the star, stay flat while fully silhouetted, and rise again as it exits.\n\n* Cosmic Ray Strikes (Red Graph): The violent, random vertical spikes in the Max Intensity plot (jumping to 80,000+) are cosmic rays—high-energy radiation from space physically striking the camera sensor.\n\n* Data Corruption (Orange & Blue Graphs): Because those cosmic rays are so blindingly bright, they mathematically drag the Mean and Standard Deviation up with them, creating artificial noise spikes in those graphs as well.\n\n* The Final Goal: The U-shape you plotted is the average of the entire image. To create the final answer key (the purple spectrum from your first image), you have to slice the image into columns and measure that exact U-shaped dip 282 separate times—once for every specific wavelength of light!","metadata":{}},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport numpy as np\n\n# Assuming clean_airs is your 3D array of shape (Time, 32, 356)\n# (If testing without the real array, we use the total length of the first axis)\ntotal_frames = clean_airs.shape[0]\nframe_indices = np.arange(total_frames)\n\n# ==========================================\n# 1. CALCULATE SUMMARY STATS PER FRAME\n# ==========================================\n# We use axis=(1, 2) to squash the Height (32) and Width (356) dimensions, \n# leaving us with exactly 1 number per time step!\n\n# Mean: The average brightness of the entire image\nmean_intensity = np.mean(clean_airs, axis=(1, 2))\n\n# Standard Deviation: The variance in pixel values (helps track contrast and background noise)\nstd_intensity = np.std(clean_airs, axis=(1, 2))\n\n# Maximum: The single brightest pixel in the entire image (perfect for finding anomalies)\nmax_intensity = np.max(clean_airs, axis=(1, 2))\n\n# ==========================================\n# 2. PLOT THE STATISTICS AS TIME SERIES\n# ==========================================\nfig, axes = plt.subplots(3, 1, figsize=(14, 10), sharex=True)\nfig.suptitle('Basic Image Statistics Over Time (AIRS-CH0)', fontsize=16, fontweight='bold')\n\n# --- Plot 1: Mean Intensity (The \"Light Curve\") ---\naxes[0].plot(frame_indices, mean_intensity, color='blue', linewidth=1)\naxes[0].set_title('1. Mean Pixel Intensity (Gradual Lighting Changes)')\naxes[0].set_ylabel('Mean Value')\naxes[0].grid(True, alpha=0.4)\n\n# --- Plot 2: Standard Deviation (Noise / Contrast) ---\naxes[1].plot(frame_indices, std_intensity, color='orange', linewidth=1)\naxes[1].set_title('2. Standard Deviation (Noise and Image Contrast)')\naxes[1].set_ylabel('Std Dev Value')\naxes[1].grid(True, alpha=0.4)\n\n# --- Plot 3: Max Intensity (Anomalies) ---\naxes[2].plot(frame_indices, max_intensity, color='red', linewidth=1)\naxes[2].set_title('3. Maximum Pixel Intensity (Sudden Spikes / Cosmic Rays)')\naxes[2].set_ylabel('Max Pixel Value')\naxes[2].set_xlabel('Frame Index (Time Step)')\naxes[2].grid(True, alpha=0.4)\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:44.891111Z","iopub.execute_input":"2026-04-25T16:37:44.891528Z","iopub.status.idle":"2026-04-25T16:37:46.887494Z","shell.execute_reply.started":"2026-04-25T16:37:44.89148Z","shell.execute_reply":"2026-04-25T16:37:46.886729Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Difference between frames (temporal dynamics)**","metadata":{}},{"cell_type":"markdown","source":"Because the \"movement\" alarm (red) and the \"shape\" alarm (green) are going off at the exact same times, it proves that the telescope slightly vibrated or shook at those exact moments.\n\nIn space terms, this is called Jitter. The spacecraft was trying to hold perfectly still for hours, but occasionally it bumped or jiggled, causing the starlight to wiggle on the camera sensor.","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport matplotlib.pyplot as plt\nimport cv2\nfrom skimage.metrics import structural_similarity as ssim\n\n# Assuming clean_airs is your 3D array of shape (Time, Height, Width)\ntotal_frames = clean_airs.shape[0]\n\npixel_diffs = []\nssim_scores = []\noptical_flow_scores = []\n\nprint(\"Calculating Temporal Dynamics... (This may take a minute for 13,000 frames!)\")\n\n# Loop through the frames and compare frame 't' to frame 't-1'\nfor t in range(1, total_frames):\n    prev_frame = clean_airs[t-1]\n    curr_frame = clean_airs[t]\n    \n    # 1. Pixel-Wise Difference (Mean Absolute Error)\n    diff = np.mean(np.abs(curr_frame - prev_frame))\n    pixel_diffs.append(diff)\n    \n    # 2. Structural Similarity Index (SSIM)\n    # We store 1 - SSIM so that higher values = bigger change\n    s_score = ssim(prev_frame, curr_frame, data_range=curr_frame.max() - curr_frame.min())\n    ssim_scores.append(1 - s_score)\n    \n    # 3. Optical Flow (Motion magnitude using OpenCV Farneback algorithm)\n    # OpenCV requires 8-bit images for standard optical flow calculations\n    prev_norm = cv2.normalize(prev_frame, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n    curr_norm = cv2.normalize(curr_frame, None, 0, 255, cv2.NORM_MINMAX).astype(np.uint8)\n    \n    flow = cv2.calcOpticalFlowFarneback(prev_norm, curr_norm, None, \n                                        pyr_scale=0.5, levels=3, winsize=15, \n                                        iterations=3, poly_n=5, poly_sigma=1.2, flags=0)\n    # Calculate magnitude of the 2D flow vectors\n    mag, ang = cv2.cartToPolar(flow[..., 0], flow[..., 1])\n    optical_flow_scores.append(np.mean(mag))\n\n# ==========================================\n# PLOT THE CHANGE SCORES OVER TIME\n# ==========================================\ntime_steps = np.arange(1, total_frames)\n\nfig, axes = plt.subplots(3, 1, figsize=(14, 10), sharex=True)\nfig.suptitle('Frame-to-Frame Temporal Dynamics', fontsize=16, fontweight='bold')\n\n# --- Plot 1: Pixel Difference ---\naxes[0].plot(time_steps, pixel_diffs, color='blue', linewidth=1)\naxes[0].set_title('1. Pixel-Wise Difference (Raw Subtraction)')\naxes[0].set_ylabel('Mean Abs Difference')\naxes[0].grid(True, alpha=0.4)\n\n# --- Plot 2: SSIM ---\naxes[1].plot(time_steps, ssim_scores, color='green', linewidth=1)\naxes[1].set_title('2. Structural Similarity (1 - SSIM)')\naxes[1].set_ylabel('Change Score')\naxes[1].grid(True, alpha=0.4)\n\n# --- Plot 3: Optical Flow ---\naxes[2].plot(time_steps, optical_flow_scores, color='red', linewidth=1)\naxes[2].set_title('3. Optical Flow Magnitude (Motion Vectors)')\naxes[2].set_ylabel('Mean Flow Vector Mag')\naxes[2].set_xlabel('Frame Index (Time Step)')\naxes[2].grid(True, alpha=0.4)\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:37:46.88854Z","iopub.execute_input":"2026-04-25T16:37:46.888847Z","iopub.status.idle":"2026-04-25T16:38:06.572891Z","shell.execute_reply.started":"2026-04-25T16:37:46.888819Z","shell.execute_reply":"2026-04-25T16:38:06.571871Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Spatial EDA (inside images)**","metadata":{}},{"cell_type":"markdown","source":"* Mostly Dark: The top histogram shows that the vast majority of pixels are nearly black.\n\n* Lighting Changes: The purple graph shows the video gets significantly darker between frames 1,000 and 4,000.\n\n* Frequent Glitches: All the tall vertical \"needles\" (spikes) across the charts indicate quick flashes or digital hiccups that happen for just a single frame.\n\n* Detail Consistency: The orange graph shows that while the video is mostly consistent, those \"glitch\" frames are much sharper/busier than the rest.\n\n* Grainy Patch: The green graph shows a dip in quality (increased noise) during that dark period in the middle of the video.","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport matplotlib.pyplot as plt\nimport cv2\nfrom scipy.ndimage import sobel\n\n# Assuming clean_airs is your 3D array of shape (Time, Height, Width)\ntotal_frames = clean_airs.shape[0]\n\n# Initialize lists to hold our spatial metrics\nmean_intensities = []\nlaplacian_vars = []\nedge_densities = []\n\nprint(\"Running Spatial EDA... (This may take a minute or two)\")\n\nfor i in range(total_frames):\n    img = clean_airs[i]\n    \n    # 1. Brightness Check (Mean Intensity)\n    mean_intensities.append(np.mean(img))\n    \n    # 2. Blur Check: Variance of Laplacian\n    # The Laplacian measures the 2nd derivative. High variance = sharp edges. Low variance = blurry.\n    lap = cv2.Laplacian(img, cv2.CV_64F)\n    laplacian_vars.append(np.var(lap))\n    \n    # 3. Sensor Noise Check: Sobel Edge Density\n    # Detects harsh edges. Extremely high values mean the image is full of sharp static.\n    sx = sobel(img, axis=0, mode='constant')\n    sy = sobel(img, axis=1, mode='constant')\n    sob_mag = np.hypot(sx, sy)\n    edge_densities.append(np.mean(sob_mag))\n\n# ==========================================\n# PLOT 1: PIXEL DISTRIBUTIONS (HISTOGRAMS)\n# ==========================================\n# We pick a few specific frames to check their internal pixel spread.\n# Feel free to change these indices to inspect different frames!\nframe_normal = 100\nframe_middle = total_frames // 2\nframe_end = total_frames - 100\n\nplt.figure(figsize=(12, 4))\nplt.hist(clean_airs[frame_normal].flatten(), bins=50, alpha=0.5, label=f'Frame {frame_normal}', color='blue')\nplt.hist(clean_airs[frame_middle].flatten(), bins=50, alpha=0.5, label=f'Frame {frame_middle}', color='green')\nplt.hist(clean_airs[frame_end].flatten(), bins=50, alpha=0.5, label=f'Frame {frame_end}', color='purple')\n\nplt.title('Pixel Distributions (Histograms) for Selected Frames', fontsize=14, fontweight='bold')\nplt.xlabel('Pixel Intensity Value')\nplt.ylabel('Number of Pixels (Frequency)')\nplt.legend()\nplt.grid(alpha=0.3)\nplt.tight_layout()\nplt.show()\n\n# ==========================================\n# PLOT 2: SPATIAL METRICS OVER TIME\n# ==========================================\ntime_steps = np.arange(total_frames)\n\nfig, axes = plt.subplots(3, 1, figsize=(14, 10), sharex=True)\nfig.suptitle('Spatial Characteristics Over Time (Inside the Images)', fontsize=16, fontweight='bold')\n\n# --- Plot A: Brightness ---\naxes[0].plot(time_steps, mean_intensities, color='purple', linewidth=1)\naxes[0].set_title('1. Mean Intensity (Checks for Underexposed/Overexposed Frames)')\naxes[0].set_ylabel('Mean Pixel Value')\naxes[0].grid(True, alpha=0.4)\n\n# --- Plot B: Sharpness/Blur ---\naxes[1].plot(time_steps, laplacian_vars, color='orange', linewidth=1)\naxes[1].set_title('2. Blur Detection (Variance of Laplacian - Drops if blurry)')\naxes[1].set_ylabel('Sharpness Score')\naxes[1].grid(True, alpha=0.4)\n\n# --- Plot C: Sensor Noise ---\naxes[2].plot(time_steps, edge_densities, color='green', linewidth=1)\naxes[2].set_title('3. Sensor Noise (Sobel Edge Density - Spikes if static/glitchy)')\naxes[2].set_ylabel('Edge Density Score')\naxes[2].set_xlabel('Frame Index (Time Step)')\naxes[2].grid(True, alpha=0.4)\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:38:06.574218Z","iopub.execute_input":"2026-04-25T16:38:06.574704Z","iopub.status.idle":"2026-04-25T16:38:10.56093Z","shell.execute_reply.started":"2026-04-25T16:38:06.574669Z","shell.execute_reply":"2026-04-25T16:38:10.559993Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Seasonality & trends (on extracted features)**","metadata":{}},{"cell_type":"markdown","source":"What These Tools Do For Your Data\n* Rolling Mean (Blue Line): The raw data is incredibly fuzzy due to quantum camera static. The rolling mean acts like a lawnmower, chopping off the static and leaving a perfectly smooth line. This makes the exoplanet's U-shaped transit completely obvious.\n\n* Rolling Variance (Orange Line): This checks for instability. If this line suddenly jumps up, it means the camera temporarily became much noisier or the spacecraft started vibrating badly.\n\nDecomposition (Green, Purple, Red): This is the ultimate un-mixing tool. It takes your complex signal and isolates it:\n\n* The Trend (Green) pulls out only the slow, massive U-shaped dip of the planet.\n\n* The Seasonality (Purple) pulls out any repeating machine rhythms (like a cooling pump turning on and off every 100 frames).\n\n* The Residuals (Red) leaves you with only the pure, unpredictable junk (like random static and sudden cosmic ray hits).","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nfrom statsmodels.tsa.seasonal import seasonal_decompose\n\n# Assuming clean_airs is your 3D array of shape (Time, Height, Width)\n# 1. Extract our 1D feature: Mean Intensity over time\nmean_intensity = np.mean(clean_airs, axis=(1, 2))\n\n# Convert to a Pandas Series for easy time-series math\nts_data = pd.Series(mean_intensity)\n\n# ==========================================\n# 2. ROLLING STATISTICS (Smoothing)\n# ==========================================\n# We use a \"window\" of 50 frames. It looks at 50 frames, averages them, \n# then moves forward by 1 frame and repeats.\nwindow_size = 50\n\n# Rolling Mean (Smooths out the random static to reveal the true shape)\nrolling_mean = ts_data.rolling(window=window_size, center=True).mean()\n\n# Rolling Variance (Measures if the \"bumpiness\" or noise is changing over time)\nrolling_var = ts_data.rolling(window=window_size, center=True).var()\n\n# ==========================================\n# 3. TIME SERIES DECOMPOSITION\n# ==========================================\n# This mathematically slices the data into 3 distinct layers:\n# - Trend (The underlying U-shape)\n# - Seasonal (Repeating rhythms/vibrations)\n# - Residual (The random noise left over)\n# We set period=100 just to look for short-term repeating patterns (like jitter).\ndecomposition = seasonal_decompose(ts_data.dropna(), model='additive', period=100)\n\ntrend = decomposition.trend\nseasonal = decomposition.seasonal\nresidual = decomposition.resid\n\n# ==========================================\n# 4. PLOTTING THE RESULTS\n# ==========================================\nfig, axes = plt.subplots(5, 1, figsize=(14, 16), sharex=True)\nfig.suptitle('Time Series Analysis: Smoothing & Decomposition', fontsize=18, fontweight='bold')\n\n# --- Plot 1: Raw Data vs. Rolling Mean ---\naxes[0].plot(ts_data, color='lightgray', label='Raw Data (Noisy)')\naxes[0].plot(rolling_mean, color='blue', linewidth=2, label=f'{window_size}-Frame Rolling Mean')\naxes[0].set_title('1. Rolling Mean (Erases static to show the true transit U-shape)')\naxes[0].set_ylabel('Intensity')\naxes[0].legend()\naxes[0].grid(True, alpha=0.4)\n\n# --- Plot 2: Rolling Variance ---\naxes[1].plot(rolling_var, color='orange', linewidth=1.5)\naxes[1].set_title('2. Rolling Variance (Spikes here mean the camera got suddenly noisy)')\naxes[1].set_ylabel('Variance')\naxes[1].grid(True, alpha=0.4)\n\n# --- Plot 3: Decomposition - TREND ---\naxes[2].plot(trend, color='green', linewidth=2)\naxes[2].set_title('3. Decomposition: TREND (The pure, slow-moving exoplanet shadow)')\naxes[2].set_ylabel('Trend')\naxes[2].grid(True, alpha=0.4)\n\n# --- Plot 4: Decomposition - SEASONALITY ---\naxes[3].plot(seasonal, color='purple', linewidth=1)\naxes[3].set_title('4. Decomposition: SEASONALITY (Repeating rhythms, e.g., thermal breathing of the telescope)')\naxes[3].set_ylabel('Seasonal')\naxes[3].grid(True, alpha=0.4)\n\n# --- Plot 5: Decomposition - RESIDUALS ---\naxes[4].plot(residual, color='red', linewidth=1)\naxes[4].set_title('5. Decomposition: RESIDUALS (Random static and unpredictable Cosmic Rays)')\naxes[4].set_ylabel('Residual')\naxes[4].set_xlabel('Frame Index (Time Step)')\naxes[4].grid(True, alpha=0.4)\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:38:10.562103Z","iopub.execute_input":"2026-04-25T16:38:10.562432Z","iopub.status.idle":"2026-04-25T16:38:12.46201Z","shell.execute_reply.started":"2026-04-25T16:38:10.562391Z","shell.execute_reply":"2026-04-25T16:38:12.461138Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Dimensionality reduction (for patterns)**","metadata":{}},{"cell_type":"markdown","source":"Summary in 3 Points\nHigh Noise: Your video is very \"grainy,\" making it hard to see changes frame-by-frame.\n\nExtreme Glitches: There are a few frames (the yellow dots) that are totally \"broken\" and don't look like the rest of the video.\n\nClear Phases: Despite the noise, the images naturally split into two distinct groups. You have successfully separated the \"normal\" parts of the video from the \"event\" parts.","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nfrom sklearn.decomposition import PCA\nfrom sklearn.manifold import TSNE\nfrom sklearn.cluster import KMeans\nimport umap  # Make sure to run: !pip install umap-learn\n\n# Assuming clean_airs is your 3D array of shape (Time, Height, Width)\nn_frames = clean_airs.shape[0]\n\nprint(\"1. Flattening images into 1D arrays...\")\n# Squash each (Height x Width) image into a single flat row of pixels\nflattened_data = clean_airs.reshape(n_frames, -1)\n\nprint(\"2. Running PCA (Fastest)...\")\npca = PCA(n_components=2)\npca_result = pca.fit_transform(flattened_data)\n\nprint(\"3. Running t-SNE (This might take a few minutes for thousands of frames!)...\")\ntsne = TSNE(n_components=2, random_state=42, init='pca', learning_rate='auto')\ntsne_result = tsne.fit_transform(flattened_data)\n\nprint(\"4. Running UMAP (Great for preserving global structure)...\")\nreducer = umap.UMAP(n_components=2, random_state=42)\numap_result = reducer.fit_transform(flattened_data)\n\nprint(\"5. Clustering similar frames (Finding Regimes)...\")\n# We use KMeans on the PCA results to group the frames into 3 distinct \"regimes\"\nkmeans = KMeans(n_clusters=3, random_state=42, n_init=10)\ncluster_labels = kmeans.fit_predict(pca_result)\n\ntime_steps = np.arange(n_frames)\n\n# ==========================================\n# PLOT A: EMBEDDING VS. TIME (Detecting Scene Changes)\n# ==========================================\n# We plot the First Component of each method over time to see the transit and outliers\nfig, axes = plt.subplots(3, 1, figsize=(14, 10), sharex=True)\nfig.suptitle('Embeddings Over Time (Component 1)', fontsize=16, fontweight='bold')\n\naxes[0].plot(time_steps, pca_result[:, 0], color='blue', linewidth=1)\naxes[0].set_title('PCA Component 1 vs Time (Watch for the U-shaped transit)')\naxes[0].grid(True, alpha=0.4)\n\naxes[1].plot(time_steps, tsne_result[:, 0], color='green', linewidth=1)\naxes[1].set_title('t-SNE Component 1 vs Time (Highly sensitive to local groupings)')\naxes[1].grid(True, alpha=0.4)\n\naxes[2].plot(time_steps, umap_result[:, 0], color='purple', linewidth=1)\naxes[2].set_title('UMAP Component 1 vs Time (Balances local and global structure)')\naxes[2].set_xlabel('Frame Index (Time)')\naxes[2].grid(True, alpha=0.4)\n\nplt.tight_layout()\nplt.show()\n\n# ==========================================\n# PLOT B: 2D CLUSTERING SCATTER PLOTS (Detecting Regimes & Outliers)\n# ==========================================\nfig, axes = plt.subplots(1, 3, figsize=(18, 5))\nfig.suptitle('2D Embedding Maps (Colored by Cluster Regime)', fontsize=16, fontweight='bold')\n\n# Helper function to plot scatters\ndef plot_scatter(ax, data, title):\n    scatter = ax.scatter(data[:, 0], data[:, 1], c=cluster_labels, cmap='viridis', s=2, alpha=0.6)\n    ax.set_title(title)\n    ax.grid(True, alpha=0.3)\n    return scatter\n\nplot_scatter(axes[0], pca_result, 'PCA Scatter')\nplot_scatter(axes[1], tsne_result, 't-SNE Scatter')\nscatter = plot_scatter(axes[2], umap_result, 'UMAP Scatter')\n\n# Add a legend for the clusters\nlegend1 = axes[2].legend(*scatter.legend_elements(), title=\"Regimes / Clusters\", loc=\"best\")\naxes[2].add_artist(legend1)\n\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:38:12.463302Z","iopub.execute_input":"2026-04-25T16:38:12.463896Z","iopub.status.idle":"2026-04-25T16:40:37.188029Z","shell.execute_reply.started":"2026-04-25T16:38:12.463849Z","shell.execute_reply":"2026-04-25T16:40:37.186912Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Detect anomalies**","metadata":{}},{"cell_type":"markdown","source":"1. Spatial Mean Intensity (The Big Event)\nThe \"Dip\": The U-shape in the middle (between frames 800 and 3800) is your actual signal—likely a transit (a planet crossing a star).\n\nThe Spikes: The vertical lines shooting up are flashes or glitches that shouldn't be there.\n\n2. Temporal Motion (The Camera Jitter)\nWhat it measures: This tracks how much the image \"shook\" or shifted from one frame to the next.\n\nWhat it tells you: Higher peaks (like the red one around frame 1200) mean the camera moved or vibrated, which can blur your data.\n\n3. Confirmed Cosmic Ray Strikes (Radiation)\nThe \"Enemy\": In space-based telescopes, high-energy particles (cosmic rays) hit the sensor and create fake bright spots.\n\nThe Red Dots: These are frames the AI is 100% sure are damaged. You usually want to delete or \"patch\" these frames because they aren't real light from the star.\n\n4. Combined ML Anomaly Score (The \"Trash\" Filter)\nThe Red Line: This is a \"danger zone\" threshold.\n\nThe Verdict: Any purple spike that crosses above that red dashed line is flagged as a bad frame. It combines the brightness, jitter, and noise to tell you: \"Don't trust this frame for science.\"","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport cv2\nfrom scipy.ndimage import sobel\nfrom statsmodels.tsa.seasonal import seasonal_decompose\nfrom sklearn.decomposition import PCA\nfrom sklearn.preprocessing import RobustScaler\nfrom sklearn.ensemble import IsolationForest\nfrom scipy.stats import zscore\nimport umap\nimport hdbscan\n\n# --- 0. PREPARATION ---\n# Assuming clean_airs is your (Time, Height, Width) array\nn = clean_airs.shape[0]\ntime_steps = np.arange(n)\n\nprint(f\"Starting Full Pipeline for {n} frames...\")\n\n# --- 1. SPATIAL EDA (Ref: image_0bdc65.png) ---\nmean_intensities = []\nlaplacian_vars = []\nedge_densities = []\n\nfor i in range(n):\n    img = clean_airs[i]\n    mean_intensities.append(np.mean(img))\n    laplacian_vars.append(np.var(cv2.Laplacian(img, cv2.CV_64F)))\n    sx = sobel(img, axis=0, mode='constant')\n    sy = sobel(img, axis=1, mode='constant')\n    edge_densities.append(np.mean(np.hypot(sx, sy)))\n\n# --- 2. TEMPORAL DYNAMICS (Ref: image_0c419c.png) ---\noptical_flow_mags = []\nfor i in range(n - 1):\n    prev, curr = clean_airs[i], clean_airs[i+1]\n    flow = cv2.calcOpticalFlowFarneback(prev, curr, None, 0.5, 3, 15, 3, 5, 1.2, 0)\n    mag, _ = cv2.cartToPolar(flow[..., 0], flow[..., 1])\n    optical_flow_mags.append(np.mean(mag))\n\n# --- 3. TIME SERIES DECOMPOSITION (Ref: image_0bca51.png) ---\nts_series = pd.Series(mean_intensities)\ndecomposition = seasonal_decompose(ts_series, model='additive', period=100)\nrolling_mean = ts_series.rolling(window=50, center=True).mean()\n\n# --- 4. DIMENSIONALITY REDUCTION (Ref: image_0b6f85.jpg) ---\nflattened = clean_airs.reshape(n, -1)\nscaled_data = RobustScaler().fit_transform(flattened)\npca_res = PCA(n_components=2).fit_transform(scaled_data)\numap_res = umap.UMAP(n_neighbors=50, min_dist=0.1).fit_transform(scaled_data)\ncluster_labels = hdbscan.HDBSCAN(min_cluster_size=30).fit_predict(umap_res)\n\n# --- 5. ANOMALY DETECTION (Ref: image_e1eff3.png) ---\nfeatures = pd.DataFrame({\n    'max_intensity': np.max(clean_airs, axis=(1, 2)),\n    'motion': [0] + list(optical_flow_mags), # Padding for length match\n    'sharpness': laplacian_vars\n})\nfeatures['z_max'] = np.abs(zscore(features['max_intensity']))\niso_forest = IsolationForest(contamination=0.015, random_state=42)\nfeatures['anomaly_label'] = iso_forest.fit_predict(features[['max_intensity', 'motion', 'sharpness']])\nanom_scores = -iso_forest.decision_function(features[['max_intensity', 'motion', 'sharpness']])\n\n# ==========================================\n# FINAL DASHBOARD GENERATION\n# ==========================================\n\n# PLOT 1: SPATIAL & TEMPORAL STATS\nfig, ax = plt.subplots(2, 1, figsize=(14, 8))\nax[0].plot(mean_intensities, color='purple', label='Mean Intensity')\nax[0].set_title('Spatial Mean Intensity (The Transit Dip)')\nax[1].plot(optical_flow_mags, color='red', label='Motion')\nax[1].set_title('Temporal Motion (Jitter)')\nplt.tight_layout()\nplt.show()\n\n# PLOT 2: THE ANOMALY ALARM (Final Report)\nfig, ax = plt.subplots(2, 1, figsize=(14, 10))\nax[0].plot(features['max_intensity'], color='gray', alpha=0.5)\ncr_idx = features.index[features['z_max'] > 5]\nax[0].scatter(cr_idx, features.loc[cr_idx, 'max_intensity'], color='red', s=20, label='Cosmic Rays')\nax[0].set_title('Confirmed Cosmic Ray Strikes')\nax[1].plot(anom_scores, color='purple')\nax[1].axhline(y=0.1, color='red', linestyle='--')\nax[1].set_title('Combined ML Anomaly Score')\nplt.show()\n\nprint(f\"Analysis Complete. Found {len(cr_idx)} cosmic rays.\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:40:37.189638Z","iopub.execute_input":"2026-04-25T16:40:37.19037Z","iopub.status.idle":"2026-04-25T16:41:32.282368Z","shell.execute_reply.started":"2026-04-25T16:40:37.190334Z","shell.execute_reply":"2026-04-25T16:41:32.281228Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Data quality checks**","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\n\ndef run_data_quality_checks(data_array):\n    print(\"🚀 Starting Data Quality Audit...\\n\")\n    \n    n_frames = data_array.shape[0]\n    report = []\n    \n    # 1. Check for Corrupted Frames (NaNs or Infs)\n    # Scientific data often has 'Not a Number' values where a sensor failed.\n    nan_count = np.isnan(data_array).sum()\n    inf_count = np.isinf(data_array).sum()\n    \n    # 2. Check Resolution Consistency\n    # We expect every frame to have the exact same Height and Width.\n    heights = [data_array[i].shape[0] for i in range(n_frames)]\n    widths = [data_array[i].shape[1] for i in range(n_frames)]\n    \n    unique_heights = set(heights)\n    unique_widths = set(widths)\n    \n    # 3. Check Format (Grayscale vs Color)\n    # For telescope data, we expect 2D (Grayscale). 3D means unexpected RGB.\n    is_grayscale = len(data_array.shape) == 3 # (Time, H, W)\n    \n    # 4. File Size / Intensity Anomalies\n    # Sudden drops to 0 mean the shutter was closed or the file didn't load.\n    frame_means = np.mean(data_array, axis=(1, 2))\n    dead_frames = np.where(frame_means <= 0)[0]\n\n    # --- SUMMARY REPORT ---\n    print(f\"--- [ QUALITY REPORT ] ---\")\n    print(f\"✅ Total Frames Scanned: {n_frames}\")\n    \n    if nan_count == 0 and inf_count == 0:\n        print(\"✅ No corrupted pixels (NaN/Inf) detected.\")\n    else:\n        print(f\"❌ WARNING: Found {nan_count} NaNs and {inf_count} Infs!\")\n\n    if len(unique_heights) == 1 and len(unique_widths) == 1:\n        print(f\"✅ Resolution is consistent: {list(unique_heights)[0]}x{list(unique_widths)[0]}\")\n    else:\n        print(f\"❌ ALERT: Resolution mismatch detected! Heights: {unique_heights}, Widths: {unique_widths}\")\n\n    if is_grayscale:\n        print(\"✅ Format Consistency: Standard Grayscale (2D per frame).\")\n    else:\n        print(f\"❌ ALERT: Format mismatch! Expected Grayscale, found {data_array.shape[3]} channels.\")\n\n    if len(dead_frames) == 0:\n        print(\"✅ No 'Dead Frames' (zero intensity) found.\")\n    else:\n        print(f\"❌ WARNING: {len(dead_frames)} frames are completely black (possible corruption).\")\n\n    # Visualize Intensity Stability\n    plt.figure(figsize=(12, 4))\n    plt.plot(frame_means, color='teal')\n    plt.title(\"Data Stability: Mean Intensity per Frame\")\n    plt.ylabel(\"Intensity\")\n    plt.xlabel(\"Frame Index\")\n    plt.grid(True, alpha=0.3)\n    plt.show()\n\n# Run the check\nrun_data_quality_checks(clean_airs)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:41:32.283635Z","iopub.execute_input":"2026-04-25T16:41:32.284559Z","iopub.status.idle":"2026-04-25T16:41:32.982593Z","shell.execute_reply.started":"2026-04-25T16:41:32.284523Z","shell.execute_reply":"2026-04-25T16:41:32.981435Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Domain-Specific Astronomy Checks**","metadata":{}},{"cell_type":"markdown","source":"1. The Graph (Star vs. Sky)\nThe Yellow Line (Star Intensity): This is the brightness of the star you are looking at. It is nice, bright, and mostly steady at the top.\n\nThe Black Line (Background): This is the empty, dark space in the corners of your image. It stays low near zero, as it should.\n\n2. The Automated Text Alerts\nThe text at the top gives you three specific diagnostic results:\n\n✅ The Good News (Star Centroid is locked): Your telescope tracking was excellent. The star stayed perfectly in focus and didn't wobble or drift off-screen.\n\n⚠️ The Warning (Background drift detected): The \"darkness\" of the empty space changed slightly (-18.23%). This usually happens when the camera sensor heats up or cools down during the recording.\n\n❌ The Alert (8 saturated pixels found): Saturated means \"overexposed.\" For 8 pixels, the star was too bright for the camera to handle, so the light maxed out the sensor. It means you lost a tiny bit of detail right at the absolute brightest center of the star.\n\nBottom line: You have very stable tracking, but your camera was slightly overexposed to the star's light, and the temperature of the sensor might have fluctuated.","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\n\ndef run_exoplanet_domain_checks(data_array):\n    print(\"🔭 Running Domain-Specific Astronomy Audit...\\n\")\n    n = data_array.shape[0]\n    \n    # 1. Background Stability Check\n    # We look at the corners of the image where there is NO star.\n    # If the background level jumps, it could be a light leak or thermal noise.\n    corners = data_array[:, :5, :5] # Top-left 5x5 patch\n    bg_stability = np.mean(corners, axis=(1, 2))\n    bg_drift = (bg_stability[-1] - bg_stability[0]) / bg_stability[0]\n    \n    # 2. Centroid Consistency (Star Drift)\n    # Does the star move around on the sensor? (Relates to image_0c419c.png)\n    # We find the \"center of mass\" for the brightest pixels.\n    drifts = []\n    for i in range(0, n, 100): # Sample every 100th frame for speed\n        img = data_array[i]\n        M = cv2.moments(img)\n        if M[\"m00\"] != 0:\n            cX = int(M[\"m10\"] / M[\"m00\"])\n            cY = int(M[\"m01\"] / M[\"m00\"])\n            drifts.append((cX, cY))\n    \n    # 3. Saturated Pixel Check\n    # If pixels hit the \"ceiling\" (e.g., 65535 for 16-bit), the data is ruined.\n    max_val = np.max(data_array)\n    sat_pixels = np.sum(data_array >= (max_val * 0.98))\n    \n    # --- DOMAIN REPORT ---\n    print(f\"--- [ ASTRONOMY DOMAIN REPORT ] ---\")\n    \n    # Background Check\n    if abs(bg_drift) < 0.05:\n        print(f\"✅ Background is stable (Drift: {bg_drift:.2%})\")\n    else:\n        print(f\"⚠️ Warning: Background drift detected ({bg_drift:.2%}). Check thermal stability.\")\n\n    # Drift Check\n    x_drift = max([d[0] for d in drifts]) - min([d[0] for d in drifts])\n    if x_drift < 2:\n        print(f\"✅ Star Centroid is locked (Drift < 2 pixels).\")\n    else:\n        print(f\"⚠️ Warning: Star is drifting by {x_drift} pixels. This will add noise to the light curve.\")\n\n    # Saturation Check\n    if sat_pixels == 0:\n        print(\"✅ No saturated pixels detected. Sensor is in linear range.\")\n    else:\n        print(f\"❌ ALERT: {sat_pixels} saturated pixels found! Data may be clipped.\")\n\n    # Visualize Background vs Star Mean\n    plt.figure(figsize=(12, 4))\n    plt.plot(bg_stability, label=\"Background (Corner)\", color='black', alpha=0.7)\n    plt.plot(np.mean(data_array, axis=(1, 2)), label=\"Star Intensity\", color='gold')\n    plt.title(\"Domain Check: Star Intensity vs Sky Background\")\n    plt.legend()\n    plt.show()\n\nrun_exoplanet_domain_checks(clean_airs)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:41:32.984074Z","iopub.execute_input":"2026-04-25T16:41:32.984522Z","iopub.status.idle":"2026-04-25T16:41:33.750158Z","shell.execute_reply.started":"2026-04-25T16:41:32.984418Z","shell.execute_reply":"2026-04-25T16:41:33.749206Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# **VIDEO**","metadata":{}},{"cell_type":"code","source":"import cv2\nimport numpy as np\n\ndef create_video_from_array(data_array, output_name='exoplanet_data.mp4', fps=30):\n    # 1. Prepare dimensions\n    n_frames, height, width = data_array.shape\n    \n    # 2. Define the codec and create VideoWriter object\n    # 'mp4v' is a standard codec that works on most players\n    fourcc = cv2.VideoWriter_fourcc(*'mp4v')\n    video = cv2.VideoWriter(output_name, fourcc, fps, (width, height), isColor=False)\n\n    print(f\"🎬 Exporting {n_frames} frames to {output_name}...\")\n\n    for i in range(n_frames):\n        # Normalize frame to 0-255 (8-bit) for video compatibility\n        frame = data_array[i]\n        frame_norm = cv2.normalize(frame, None, 0, 255, cv2.NORM_MINMAX).astype('uint8')\n        \n        video.write(frame_norm)\n\n    video.release()\n    print(\"✅ Video saved successfully!\")\n\n# Run it on your dataset\ncreate_video_from_array(clean_airs)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:41:44.929628Z","iopub.execute_input":"2026-04-25T16:41:44.930088Z","iopub.status.idle":"2026-04-25T16:41:45.717381Z","shell.execute_reply.started":"2026-04-25T16:41:44.930045Z","shell.execute_reply":"2026-04-25T16:41:45.716436Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nfrom matplotlib.animation import FFMpegWriter\n\n# 1. Setup the figure\nfig, (ax1, ax2) = plt.subplots(1, 2, figsize=(15, 6))\nimg_display = ax1.imshow(clean_airs[0], cmap='viridis')\nax1.set_title(\"Telescope View\")\n\n# Light curve setup\nline, = ax2.plot([], [], color='purple', lw=2)\nax2.set_xlim(0, len(mean_intensities))\nax2.set_ylim(min(mean_intensities), max(mean_intensities))\nax2.set_title(\"Real-time Light Curve\")\nax2.set_xlabel(\"Frame Index\")\n\n# 2. Update function for each frame\ndef update(i):\n    img_display.set_data(clean_airs[i])\n    line.set_data(range(i), mean_intensities[:i])\n    return img_display, line\n\n# 3. Save the animation\nmetadata = dict(title='Exoplanet Analysis', artist='Matplotlib')\nwriter = FFMpegWriter(fps=20, metadata=metadata)\n\nwith writer.saving(fig, \"analysis_dashboard.mp4\", dpi=100):\n    for i in range(0, len(clean_airs), 5): # Skip by 5 for faster encoding\n        update(i)\n        writer.grab_frame()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:41:50.696938Z","iopub.execute_input":"2026-04-25T16:41:50.697397Z","iopub.status.idle":"2026-04-25T16:43:06.112311Z","shell.execute_reply.started":"2026-04-25T16:41:50.697354Z","shell.execute_reply":"2026-04-25T16:43:06.111432Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nfrom matplotlib.animation import FuncAnimation\nfrom IPython.display import HTML\n\n# 1. Setup the figure and subplots\nfig, (ax_img, ax_plot) = plt.subplots(1, 2, figsize=(16, 6))\nfig.suptitle('Synchronized Exoplanet Transit Analysis', fontsize=16, fontweight='bold')\n\n# Left Side: The Image\nimg_display = ax_img.imshow(clean_airs[0], cmap='viridis', aspect='auto')\nax_img.set_title(\"Telescope Raw Frame\")\nplt.colorbar(img_display, ax=ax_img, label='Intensity')\n\n# Right Side: The Light Curve\nax_plot.plot(mean_intensities, color='gray', alpha=0.3, label='Full Light Curve')\nax_plot.set_title(\"Intensity Over Time (The Transit)\")\nax_plot.set_xlabel(\"Frame Index\")\nax_plot.set_ylabel(\"Mean Intensity\")\n\n# The \"Moving Cursor\" - a vertical line that follows the current frame\ncursor = ax_plot.axvline(x=0, color='red', linestyle='--', linewidth=2, label='Current Frame')\n# A small dot that follows the exact data point\npoint, = ax_plot.plot([], [], 'ro', markersize=8) \nax_plot.legend()\n\n# 2. The Animation Function\ndef update(i):\n    # Update the image\n    img_display.set_data(clean_airs[i])\n    \n    # Update the cursor position on the plot\n    cursor.set_xdata([i])\n    \n    # Update the red dot on the plot\n    point.set_data([i], [mean_intensities[i]])\n    \n    return img_display, cursor, point\n\n# 3. Create the Animation\n# 'frames' is every 10th frame for speed; change to len(clean_airs) for full detail\nani = FuncAnimation(fig, update, frames=range(0, len(clean_airs), 10), interval=50, blit=True)\n\n# 4. SAVE AND DISPLAY\n# This saves the file to your current folder\nani.save('synchronized_transit_dashboard.mp4', writer='ffmpeg')\n\n# This shows the video RIGHT HERE in your notebook cell\nplt.close() # Prevents extra static plot from showing\nHTML(ani.to_html5_video())","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:45:40.165043Z","iopub.execute_input":"2026-04-25T16:45:40.165447Z","iopub.status.idle":"2026-04-25T16:49:52.828384Z","shell.execute_reply.started":"2026-04-25T16:45:40.16541Z","shell.execute_reply":"2026-04-25T16:49:52.827379Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nfrom matplotlib.animation import FuncAnimation\nfrom matplotlib.gridspec import GridSpec\nfrom IPython.display import HTML\nimport numpy as np\nimport pandas as pd\n\n# --- PREP DATA (Assuming features DataFrame exists from step 9) ---\nn = len(clean_airs)\ntime_indices = np.arange(n)\n# Calculate a very smooth rolling mean for the main presentation plot\nsmoothed_curve = pd.Series(mean_intensities).rolling(window=100, center=True).mean().fillna(method='bfill').fillna(method='ffill')\nanomaly_scores = anom_scores # From Isolation Forest decission_function\n\n# ==========================================\n# 1. SETUP THE DASHBOARD LAYOUT\n# ==========================================\nfig = plt.figure(figsize=(18, 10))\ngs = GridSpec(2, 2, figure=fig, width_ratios=[1, 1.2], height_ratios=[1, 1])\n\n# A. Raw Telescope Image\nax_img = fig.add_subplot(gs[:, 0])\nimg_display = ax_img.imshow(clean_airs[0], cmap='viridis', aspect='auto', interpolation='nearest')\nax_img.set_title(\"Telescope View (Raw Data)\", fontsize=14, fontweight='bold')\nplt.colorbar(img_display, ax=ax_img, label='Pixel Intensity')\n\n# B. Polished Light Curve (The Science)\nax_lc = fig.add_subplot(gs[0, 1])\nax_lc.plot(time_indices, smoothed_curve, color='teal', linewidth=3, label='Smoothed Flux')\nax_lc.set_title(\"Scientific Signal (Normalized Transit Curve)\", fontsize=14, fontweight='bold')\nax_lc.set_ylabel(\"Relative Brightness\")\nax_lc.grid(True, alpha=0.3)\n# Moving cursor on LC\ncursor_lc = ax_lc.axvline(x=0, color='red', linestyle='--', linewidth=2)\nlabel_ann = ax_lc.annotate('', xy=(0, min(smoothed_curve)), xytext=(20, 20),\n                           textcoords='offset points', bbox=dict(boxstyle='round', fc='w', alpha=0.8),\n                           arrowprops=dict(arrowstyle='->'))\n\n# C. Anomaly Alarms (The Math)\nax_anom = fig.add_subplot(gs[1, 1], sharex=ax_lc)\nax_anom.fill_between(time_indices, -anomaly_scores, 0.1, color='purple', alpha=0.2, label='Normal Range')\nax_anom.fill_between(time_indices, -anomaly_scores, 0.1, where=(-anomaly_scores > 0.1), color='red', alpha=0.6, label='ALERT (Cosmic Ray/Jitter)')\nax_anom.set_title(\"Automated Data Quality Monitor\", fontsize=14, fontweight='bold')\nax_anom.set_xlabel(\"Frame Index (Time)\")\nax_anom.set_ylabel(\"Anomaly Score\")\nax_anom.legend(loc='upper right', fontsize=8)\n# Moving cursor on Anomaly Plot\ncursor_anom = ax_anom.axvline(x=0, color='red', linestyle='--', linewidth=2)\n\nfig.tight_layout(rect=[0, 0.03, 1, 0.95])\n\n# ==========================================\n# 2. THE ANIMATION FUNCTION\n# ==========================================\n# Skip frames for faster rendering (Render every 20th frame)\nframes_to_render = range(0, n, 20)\n\ndef update_dashboard(i):\n    # Update Image\n    img_display.set_data(clean_airs[i])\n    # Optionally dynamically adjust color scaling to highlight flicker\n    # img_display.set_clim(np.percentile(clean_airs[i], 1), np.percentile(clean_airs[i], 99))\n    \n    # Update Cursors\n    cursor_lc.set_xdata([i])\n    cursor_anom.set_xdata([i])\n    \n    # Update Dynamic Label (Example labels based on typical frame indices)\n    if i < 1500: status = \"Baseline (Star)\"\n    elif 1500 <= i < 1800: status = \"Ingress (Transit Starts)\"\n    elif 1800 <= i < 3500: status = \"In-Transit (Planet Shadow)\"\n    elif 3500 <= i < 3800: status = \"Egress (Transit Ends)\"\n    else: status = \"Baseline\"\n    \n    label_ann.set_text(f\"Frame {i}: {status}\")\n    label_ann.xy = (i, smoothed_curve[i])\n\n    return img_display, cursor_lc, cursor_anom, label_ann\n\n# ==========================================\n# 3. GENERATE AND SAVE\n# ==========================================\nprint(\"🎬 Rendering presentation dashboard...\")\nani_dash = FuncAnimation(fig, update_dashboard, frames=frames_to_render, blit=True)\n\n# Save as polished MP4\nani_dash.save('exoplanet_executive_dashboard.mp4', writer='ffmpeg', fps=30, dpi=150)\nplt.close()\n\nprint(\"✅ Complete. Displaying preview below:\")\n# Display directly in notebook for checking\nHTML(ani_dash.to_html5_video())","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-04-25T16:51:55.510916Z","iopub.execute_input":"2026-04-25T16:51:55.511436Z","iopub.status.idle":"2026-04-25T16:58:09.393088Z","shell.execute_reply.started":"2026-04-25T16:51:55.511391Z","shell.execute_reply":"2026-04-25T16:58:09.391695Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}