{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.11","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":39763,"databundleVersionId":11756775,"sourceType":"competition"}],"dockerImageVersionId":31012,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import numpy as np\nimport tensorflow as tf\nfrom tensorflow.keras import layers, models, optimizers\nfrom pathlib import Path\nimport matplotlib.pyplot as plt\nfrom sklearn.model_selection import train_test_split\nfrom tqdm import tqdm\nimport os\n\nnp.random.seed(42)\ntf.random.set_seed(42)\n\nBATCH_SIZE = 4\nEPOCHS = 30\nLEARNING_RATE = 1e-4\nWAVELET_SCALE = 1e-3\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-04-22T21:50:48.597672Z","iopub.execute_input":"2025-04-22T21:50:48.598598Z","iopub.status.idle":"2025-04-22T21:50:48.632358Z","shell.execute_reply.started":"2025-04-22T21:50:48.598566Z","shell.execute_reply":"2025-04-22T21:50:48.631397Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def get_data_paths(data_dir):\n    input_files = []\n    output_files = []\n    for root, _, files in os.walk(data_dir):\n        for file in files:\n            if file.endswith('.npy'):\n                full_path = os.path.join(root, file)\n                if 'data' in file or 'seis' in file:\n                    input_files.append(full_path)\n                    output_file = full_path.replace('data', 'model').replace('seis', 'vel')\n                    if os.path.exists(output_file):\n                        output_files.append(output_file)\n    return input_files, output_files\n\ndef load_and_preprocess_data(input_files, output_files):\n    X, y = [], []\n    for inp, out in tqdm(zip(input_files, output_files), total=len(input_files)):\n        try:\n            seismic = np.load(inp)\n            velocity = np.load(out)\n            seismic = np.mean(seismic[0, :, :, :], axis=0)\n            velocity = velocity[0, 0, :, :]\n            seismic = (seismic - np.mean(seismic)) / (np.std(seismic) + 1e-8) * WAVELET_SCALE\n            velocity = (velocity - 1500) / 4500\n            X.append(seismic)\n            y.append(velocity)\n        except Exception as e:\n            print(f\"Error loading {inp} or {out}: {str(e)}\")\n    return np.array(X), np.array(y)\n\ntrain_data_dir = '/kaggle/input/waveform-inversion/train_samples'\ntest_data_dir = '/kaggle/input/waveform-inversion/test'\n\nall_inputs, all_outputs = get_data_paths(train_data_dir)\n\nX, y = load_and_preprocess_data(all_inputs, all_outputs)\n\nX_train, X_val, y_train, y_val = train_test_split(X, y, test_size=0.2, random_state=42)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-04-22T21:50:48.634696Z","iopub.execute_input":"2025-04-22T21:50:48.635020Z","iopub.status.idle":"2025-04-22T21:50:54.086911Z","shell.execute_reply.started":"2025-04-22T21:50:48.634999Z","shell.execute_reply":"2025-04-22T21:50:54.085899Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def build_pinn_model(input_shape, output_shape):\n    inputs = layers.Input(shape=input_shape)\n    x = layers.Conv1D(32, 3, activation='relu', padding='same')(inputs)\n    x = layers.MaxPooling1D(2)(x)\n    x = layers.Conv1D(64, 3, activation='relu', padding='same')(x)\n    x = layers.MaxPooling1D(2)(x)\n    x = layers.Flatten()(x)\n    x = layers.Dense(256, activation='relu')(x)\n    x = layers.Dense(512, activation='relu')(x)\n    x = layers.Dense(np.prod(output_shape), activation='relu')(x)\n    outputs = layers.Reshape(output_shape)(x)\n    model = models.Model(inputs, outputs)\n    return model\n\ninput_shape = X_train.shape[1:]\noutput_shape = y_train.shape[1:]\nmodel = build_pinn_model(input_shape, output_shape)\nmodel.summary()\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-04-22T21:50:54.087877Z","iopub.execute_input":"2025-04-22T21:50:54.088158Z","iopub.status.idle":"2025-04-22T21:50:54.228311Z","shell.execute_reply.started":"2025-04-22T21:50:54.088136Z","shell.execute_reply":"2025-04-22T21:50:54.227539Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def physics_loss(y_true, y_pred):\n    data_loss = tf.reduce_mean(tf.abs(y_true - y_pred))\n    y_pred_4d = tf.expand_dims(y_pred, axis=-1)\n    dy_dx = tf.image.image_gradients(y_pred_4d)\n    physics_constraint = tf.reduce_mean(tf.abs(dy_dx[0]) + tf.reduce_mean(tf.abs(dy_dx[1])))\n    return data_loss + 0.1 * physics_constraint\n\nmodel.compile(optimizer=optimizers.Adam(LEARNING_RATE),\n              loss=physics_loss,\n              metrics=['mae'])\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-04-22T21:50:54.230038Z","iopub.execute_input":"2025-04-22T21:50:54.230305Z","iopub.status.idle":"2025-04-22T21:50:54.248505Z","shell.execute_reply.started":"2025-04-22T21:50:54.230286Z","shell.execute_reply":"2025-04-22T21:50:54.247763Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"callbacks = [\n    tf.keras.callbacks.ReduceLROnPlateau(patience=5),\n    tf.keras.callbacks.EarlyStopping(patience=10, restore_best_weights=True),\n    tf.keras.callbacks.ModelCheckpoint('best_model.keras', save_best_only=True)\n]\n\nhistory = model.fit(\n    X_train, y_train,\n    validation_data=(X_val, y_val),\n    batch_size=BATCH_SIZE,\n    epochs=EPOCHS,\n    callbacks=callbacks\n)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-04-22T21:50:54.249311Z","iopub.execute_input":"2025-04-22T21:50:54.249570Z","iopub.status.idle":"2025-04-22T21:51:11.571510Z","shell.execute_reply.started":"2025-04-22T21:50:54.249546Z","shell.execute_reply":"2025-04-22T21:51:11.570581Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"plt.figure(figsize=(12, 4))\nplt.subplot(1, 2, 1)\nplt.plot(history.history['loss'], label='train')\nplt.plot(history.history['val_loss'], label='val')\nplt.legend()\nplt.title('Loss')\n\nplt.subplot(1, 2, 2)\nplt.plot(history.history['mae'], label='train')\nplt.plot(history.history['val_mae'], label='val')\nplt.legend()\nplt.title('MAE')\nplt.show()\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-04-22T21:51:11.573259Z","iopub.execute_input":"2025-04-22T21:51:11.573831Z","iopub.status.idle":"2025-04-22T21:51:11.951201Z","shell.execute_reply.started":"2025-04-22T21:51:11.573796Z","shell.execute_reply":"2025-04-22T21:51:11.950184Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def generate_submission(model, test_dir):\n    test_files = [f for f in Path(test_dir).rglob('*.npy')]\n    x_cols = [f\"x_{i}\" for i in range(1, 70, 2)]\n    header = [\"oid_ypos\"] + x_cols\n    submission = [\",\".join(header)]\n    \n    for file in tqdm(test_files):\n        oid = file.stem\n        data = np.load(file)\n\n        if data.ndim == 4:\n            data = data[0]\n        if data.ndim == 3:\n            data = np.mean(data, axis=0)\n\n        data = (data - np.mean(data)) / (np.std(data) + 1e-8) * WAVELET_SCALE\n\n        pred = model.predict(np.expand_dims(data, 0))[0]\n        pred = pred * 4500 + 1500  \n\n        for y_pos in range(70):\n            odd_cols = pred[y_pos, 1::2]\n\n            line = f\"{oid}_y_{y_pos},\" + \",\".join([f\"{v:.1f}\" for v in odd_cols])\n            submission.append(line)\n\n    with open('submission.csv', 'w') as f:\n        f.write(\"\\n\".join(submission))\n\ngenerate_submission(model, test_data_dir)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-04-22T21:52:58.636801Z","iopub.execute_input":"2025-04-22T21:52:58.637282Z"}},"outputs":[],"execution_count":null}]}