Compare commits
3 Commits
dev_01
..
a4ad11742a
| Author | SHA1 | Date | |
|---|---|---|---|
| a4ad11742a | |||
| f2c11079e4 | |||
| b5b600023e |
+6
@@ -0,0 +1,6 @@
|
|||||||
|
<component name="InspectionProjectProfileManager">
|
||||||
|
<settings>
|
||||||
|
<option name="USE_PROJECT_PROFILE" value="false" />
|
||||||
|
<version value="1.0" />
|
||||||
|
</settings>
|
||||||
|
</component>
|
||||||
Generated
+8
@@ -0,0 +1,8 @@
|
|||||||
|
<?xml version="1.0" encoding="UTF-8"?>
|
||||||
|
<project version="4">
|
||||||
|
<component name="ProjectModuleManager">
|
||||||
|
<modules>
|
||||||
|
<module fileurl="file://$PROJECT_DIR$/.idea/remote-sensing.iml" filepath="$PROJECT_DIR$/.idea/remote-sensing.iml" />
|
||||||
|
</modules>
|
||||||
|
</component>
|
||||||
|
</project>
|
||||||
Generated
+8
@@ -0,0 +1,8 @@
|
|||||||
|
<?xml version="1.0" encoding="UTF-8"?>
|
||||||
|
<module type="PYTHON_MODULE" version="4">
|
||||||
|
<component name="NewModuleRootManager">
|
||||||
|
<content url="file://$MODULE_DIR$" />
|
||||||
|
<orderEntry type="inheritedJdk" />
|
||||||
|
<orderEntry type="sourceFolder" forTests="false" />
|
||||||
|
</component>
|
||||||
|
</module>
|
||||||
Generated
+6
@@ -0,0 +1,6 @@
|
|||||||
|
<?xml version="1.0" encoding="UTF-8"?>
|
||||||
|
<project version="4">
|
||||||
|
<component name="VcsDirectoryMappings">
|
||||||
|
<mapping directory="" vcs="Git" />
|
||||||
|
</component>
|
||||||
|
</project>
|
||||||
Generated
+48
@@ -0,0 +1,48 @@
|
|||||||
|
<?xml version="1.0" encoding="UTF-8"?>
|
||||||
|
<project version="4">
|
||||||
|
<component name="ChangeListManager">
|
||||||
|
<list default="true" id="8512e5bb-2b73-4d09-a37b-d0b357b3fbe9" name="Changes" comment="" />
|
||||||
|
<option name="SHOW_DIALOG" value="false" />
|
||||||
|
<option name="HIGHLIGHT_CONFLICTS" value="true" />
|
||||||
|
<option name="HIGHLIGHT_NON_ACTIVE_CHANGELIST" value="false" />
|
||||||
|
<option name="LAST_RESOLUTION" value="IGNORE" />
|
||||||
|
</component>
|
||||||
|
<component name="Git.Settings">
|
||||||
|
<option name="RECENT_GIT_ROOT_PATH" value="$PROJECT_DIR$" />
|
||||||
|
</component>
|
||||||
|
<component name="ProjectColorInfo"><![CDATA[{
|
||||||
|
"associatedIndex": 1
|
||||||
|
}]]></component>
|
||||||
|
<component name="ProjectId" id="39nRCQRaBb6bqtrPoIoBjpMe8Bs" />
|
||||||
|
<component name="ProjectViewState">
|
||||||
|
<option name="hideEmptyMiddlePackages" value="true" />
|
||||||
|
<option name="showLibraryContents" value="true" />
|
||||||
|
</component>
|
||||||
|
<component name="PropertiesComponent"><![CDATA[{
|
||||||
|
"keyToString": {
|
||||||
|
"ModuleVcsDetector.initialDetectionPerformed": "true",
|
||||||
|
"RunOnceActivity.ShowReadmeOnStart": "true",
|
||||||
|
"RunOnceActivity.TerminalTabsStorage.copyFrom.TerminalArrangementManager.252": "true",
|
||||||
|
"RunOnceActivity.git.unshallow": "true",
|
||||||
|
"git-widget-placeholder": "dev__01",
|
||||||
|
"last_opened_file_path": "//wsl.localhost/Ubuntu-22.04/home/x79/remote-sensing"
|
||||||
|
}
|
||||||
|
}]]></component>
|
||||||
|
<component name="SharedIndexes">
|
||||||
|
<attachedChunks>
|
||||||
|
<set>
|
||||||
|
<option value="bundled-python-sdk-4762d8aabb82-6d6dccd035ac-com.jetbrains.pycharm.pro.sharedIndexes.bundled-PY-253.30387.173" />
|
||||||
|
</set>
|
||||||
|
</attachedChunks>
|
||||||
|
</component>
|
||||||
|
<component name="TaskManager">
|
||||||
|
<task active="true" id="Default" summary="Default task">
|
||||||
|
<changelist id="8512e5bb-2b73-4d09-a37b-d0b357b3fbe9" name="Changes" comment="" />
|
||||||
|
<created>1771329720024</created>
|
||||||
|
<option name="number" value="Default" />
|
||||||
|
<option name="presentableId" value="Default" />
|
||||||
|
<updated>1771329720024</updated>
|
||||||
|
</task>
|
||||||
|
<servers />
|
||||||
|
</component>
|
||||||
|
</project>
|
||||||
File diff suppressed because one or more lines are too long
File diff suppressed because one or more lines are too long
File diff suppressed because one or more lines are too long
@@ -1,171 +0,0 @@
|
|||||||
#!/usr/bin/env python
|
|
||||||
# coding: utf-8
|
|
||||||
|
|
||||||
# In[1]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '%matplotlib inline\n\nimport importlib\nimport new_import_ODC \n\nimportlib.reload(new_import_ODC)\n\nfrom new_import_ODC import *\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[2]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Cấu hình Daskgateway\ncluster, client = notebook_utils.initialize_dask(use_gateway=True, workers=(1, 10))\n# Khai báo 1 Datacube là dc\ndc = None\n\n# Cấu hình truy cập dịch vụ S3\nconfigure_s3_access(aws_unsigned=False, requester_pays=True, client=client)\n\nclient\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[3]:
|
|
||||||
|
|
||||||
|
|
||||||
## cấu hình thời gian lấy ảnh và tọa độ
|
|
||||||
date_range = ("2022-09-01", "2022-10-01")
|
|
||||||
longtitude_range = (105.86, 105.94)
|
|
||||||
latitude_range = (9.65, 9.69)
|
|
||||||
|
|
||||||
coordinates = (longtitude_range, latitude_range)
|
|
||||||
|
|
||||||
|
|
||||||
# In[4]:
|
|
||||||
|
|
||||||
|
|
||||||
## truy vấn ảnh vệ tinh sen2
|
|
||||||
data = load_data(None, date_range, longtitude_range, latitude_range)
|
|
||||||
notebook_utils.heading(notebook_utils.xarray_object_size(data))
|
|
||||||
display(data)
|
|
||||||
|
|
||||||
|
|
||||||
# In[5]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Tiến hành loại bỏ các vị trí bị mây ảnh hưởng\nresult = mask_clean(data)\n# progress(result)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[6]:
|
|
||||||
|
|
||||||
|
|
||||||
# Tiến hành tính toán NDVI
|
|
||||||
ds1 = calculate_indices(result, index="NDVI", satellite_mission="s2")
|
|
||||||
ndvi = ds1["NDVI"]
|
|
||||||
display(ndvi)
|
|
||||||
|
|
||||||
|
|
||||||
# In[7]:
|
|
||||||
|
|
||||||
|
|
||||||
## Hiển thị ảnh NDVI chưa điền các giá trị mây (chưa fill nan)
|
|
||||||
plt.imshow(ndvi.isel(time=0))
|
|
||||||
|
|
||||||
|
|
||||||
# In[8]:
|
|
||||||
|
|
||||||
|
|
||||||
# Thiết lập giá trị trung bình mùa vụ để xử lý các điểm ảnh bị mây dựa vào sự thay đổi theo mùa
|
|
||||||
time_split = [
|
|
||||||
slice("2022-09-01", "2023-01-01"),
|
|
||||||
slice("2023-01-01", "2023-05-01"),
|
|
||||||
slice("2023-05-01", "2023-07-01"),
|
|
||||||
slice("2023-07-01", "2022-10-01"),
|
|
||||||
]
|
|
||||||
|
|
||||||
# Điền mây ở các vị trí mang giá trị nan (fill nan)
|
|
||||||
fill_nan_ndvi = fill_nan(ndvi, time_split)
|
|
||||||
|
|
||||||
# In kết quả ảnh NDVI đã điền mây (đã fill nan)
|
|
||||||
plt.imshow(fill_nan_ndvi.isel(time=0))
|
|
||||||
|
|
||||||
|
|
||||||
# In[9]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '## tính ndvi theo tháng\naverage_ndvi = fill_nan_ndvi.resample(time="1M").mean().persist()\n# progress(average_ndvi)\n\n# compute average_ndvi\naverage_ndvi = average_ndvi.compute()\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[10]:
|
|
||||||
|
|
||||||
|
|
||||||
#Load dữ liệu ảnh Sentinel 1
|
|
||||||
dsvh, dsvv = load_data_sen1(None, date_range, coordinates)
|
|
||||||
average_vv = calculate_average(dsvv, time_pattern='1M')
|
|
||||||
average_vh = calculate_average(dsvh, time_pattern='1M')
|
|
||||||
|
|
||||||
|
|
||||||
# In[11]:
|
|
||||||
|
|
||||||
|
|
||||||
## cấu hình bộ dữ liệu điểm huấn luyện mô hình (train file)
|
|
||||||
train_path = "train/ST_training_data_updated_1130points_new.shp" # đường dẫn shp file train
|
|
||||||
|
|
||||||
## load dữ liệu điểm huấn luyện mô hình (train file)
|
|
||||||
train = load_train_data(train_path)
|
|
||||||
train.head()
|
|
||||||
|
|
||||||
# cấu hình nhãn dữ liệu
|
|
||||||
label_mapping = {
|
|
||||||
"Lua tom": "0",
|
|
||||||
"Lua": "1",
|
|
||||||
"CHN": "2",
|
|
||||||
"CLN": "3",
|
|
||||||
"TS": "4",
|
|
||||||
"Song": "5",
|
|
||||||
"Dat xay dung": "6",
|
|
||||||
"Rung": "7",
|
|
||||||
}
|
|
||||||
|
|
||||||
# xây dựng tập dữ liệu (dataset) chứa dữ liệu VH, VV, NDVI
|
|
||||||
datasets = get_data_sen1_and_sen2(train, average_ndvi, average_vh, average_vv)
|
|
||||||
|
|
||||||
# chia tập dữ liệu thành các phần theo tỉ lệ 80(80-20)-20 tương ứng với tập train, validate, test
|
|
||||||
X_train, X_val, X_test, y_train, y_val, y_test = split_train_data(
|
|
||||||
train, label_mapping, datasets
|
|
||||||
)
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Import XGBoost\nimport xgboost as xgb\nfrom sklearn.metrics import accuracy_score\nimport numpy as np\n\n# Convert to numpy arrays\nX_train_np = np.asarray(X_train, dtype=np.float32)\nX_val_np = np.asarray(X_val, dtype=np.float32)\ny_train_np = np.asarray(y_train, dtype=np.int32)\ny_val_np = np.asarray(y_val, dtype=np.int32)\n\nprint("🚀 Training XGBoost model...")\nprint(f" Train samples: {len(X_train_np)}")\nprint(f" Val samples: {len(X_val_np)}")\nprint(f" Features: {X_train_np.shape[1]}")\nprint(f" Classes: 8\\n")\n\n# XGBoost parameters\nparams = {\n \'objective\': \'multi:softmax\', # Multi-class classification\n \'num_class\': 8, # 8 land use classes\n \'max_depth\': 6, # Maximum tree depth\n \'learning_rate\': 0.1, # Learning rate\n \'n_estimators\': 200, # Number of trees\n \'subsample\': 0.8, # Subsample ratio\n \'colsample_bytree\': 0.8, # Feature sampling ratio\n \'random_state\': 42,\n \'n_jobs\': -1, # Use all CPU cores\n \'eval_metric\': \'mlogloss\' # Multi-class log loss\n}\n\n# Train XGBoost model\nmodel = xgb.XGBClassifier(**params)\n\nmodel.fit(\n X_train_np, y_train_np,\n eval_set=[(X_train_np, y_train_np), (X_val_np, y_val_np)],\n verbose=True\n)\n\n# Validation accuracy\ny_val_pred = model.predict(X_val_np)\nval_accuracy = accuracy_score(y_val_np, y_val_pred)\nprint(f"\\n✅ Training completed!")\nprint(f" Validation Accuracy: {val_accuracy:.4f} ({val_accuracy*100:.2f}%)")\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Evaluate on test set\nX_test_np = np.asarray(X_test, dtype=np.float32)\ny_test_np = np.asarray(y_test, dtype=np.int32)\n\nprint("📊 Evaluating XGBoost model on test set...\\n")\n\n# Predictions\ny_pred_test = model.predict(X_test_np)\n\n# Metrics\nfrom sklearn.metrics import accuracy_score, precision_score, recall_score, f1_score, confusion_matrix\n\ntest_accuracy = accuracy_score(y_test_np, y_pred_test)\nprecision = precision_score(y_test_np, y_pred_test, average=\'weighted\', zero_division=0)\nrecall = recall_score(y_test_np, y_pred_test, average=\'weighted\', zero_division=0)\nf1 = f1_score(y_test_np, y_pred_test, average=\'weighted\', zero_division=0)\n\nprint(f"📈 Test Results:")\nprint(f" Accuracy: {test_accuracy:.4f} ({test_accuracy*100:.2f}%)")\nprint(f" Precision: {precision:.4f}")\nprint(f" Recall: {recall:.4f}")\nprint(f" F1-Score: {f1:.4f}\\n")\n\n# Confusion Matrix\nfrom sklearn.metrics import ConfusionMatrixDisplay\nimport matplotlib.pyplot as plt\n\n# Create figure first\nfig, ax = plt.subplots(figsize=(10, 8))\n\nclass_names = list(label_mapping.keys())\ncm = confusion_matrix(y_test_np, y_pred_test)\ndisp = ConfusionMatrixDisplay(confusion_matrix=cm, display_labels=class_names)\ndisp.plot(cmap=\'Blues\', ax=ax)\nplt.xticks(rotation=45, ha=\'right\')\nplt.title(\'XGBoost Confusion Matrix\')\nplt.tight_layout()\nplt.show()\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
# Lưu mô hình huấn luyện
|
|
||||||
import json
|
|
||||||
import joblib
|
|
||||||
|
|
||||||
# Save XGBoost model
|
|
||||||
model_path = "model_xgboost.joblib"
|
|
||||||
joblib.dump(model, model_path)
|
|
||||||
print(f"✅ Model saved to {model_path}")
|
|
||||||
|
|
||||||
# Save model info
|
|
||||||
info = {
|
|
||||||
"model_type": "XGBoost",
|
|
||||||
"num_classes": 8,
|
|
||||||
"classes": list(label_mapping.keys()),
|
|
||||||
"num_features": X_train_np.shape[1],
|
|
||||||
"params": params,
|
|
||||||
"accuracy": float(test_accuracy),
|
|
||||||
"precision": float(precision),
|
|
||||||
"recall": float(recall),
|
|
||||||
"f1_score": float(f1),
|
|
||||||
}
|
|
||||||
|
|
||||||
with open("model_xgboost_info.json", "w") as f:
|
|
||||||
json.dump(info, f, indent=2)
|
|
||||||
|
|
||||||
print(f"✅ Model info saved to model_xgboost_info.json")
|
|
||||||
|
|
||||||
|
|
||||||
# In[15]:
|
|
||||||
|
|
||||||
|
|
||||||
# đóng client, cluster
|
|
||||||
# client.close()
|
|
||||||
# cluster.close()
|
|
||||||
|
|
||||||
@@ -1,56 +0,0 @@
|
|||||||
#!/usr/bin/env python
|
|
||||||
# coding: utf-8
|
|
||||||
|
|
||||||
# In[6]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '%matplotlib inline\n\n# Import Microsoft Planetary Computer libraries\nimport planetary_computer\nfrom pystac_client import Client\nfrom odc.stac import load as stac_load\n\n# Standard imports\nimport xarray as xr\nimport numpy as np\nimport matplotlib.pyplot as plt\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.metrics import accuracy_score, classification_report, confusion_matrix, ConfusionMatrixDisplay\nimport geopandas as gpd\n\n# XGBoost for GPU training\nimport xgboost as xgb\n\nfrom xgboost import XGBClassifier\n\nprint(f" XGBoost version: {xgb.__version__}")\n\nprint("✅ All modules loaded successfully")\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[7]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Kết nối tới Microsoft Planetary Computer STAC\nfrom pystac_client import Client\n\n# KHÔNG dùng modifier ở catalog level để tránh items bị convert thành dict\ncatalog = Client.open(\n "https://planetarycomputer.microsoft.com/api/stac/v1"\n)\nprint("✅ Connected to Microsoft Planetary Computer")\n\nprint("\\n" + "="*70)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[8]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# 🌍 Định nghĩa khu vực và thời gian\nprint("="*70)\nprint("CONFIGURATION")\nprint("="*70)\n\n# Khu vực quan tâm (Vietnam - Mekong Delta) - GIẢM DIỆN TÍCH ~40%\nbbox = [105.6, 9.3, 106.2, 9.8] # [min_lon, min_lat, max_lon, max_lat]\n\n# GIẢM THỜI GIAN xuống 3 tháng để giảm kích thước dữ liệu cho PC\ntime_range = "2023-03-01/2023-05-31" # 3 tháng (mùa khô)\n\nprint(f"\\n📍 Area of Interest:")\nprint(f" Longitude: {bbox[0]} to {bbox[2]}")\nprint(f" Latitude: {bbox[1]} to {bbox[3]}")\nprint(f"\\n📅 Time Range: {time_range}")\nprint(f" ⚠️ Optimized for personal computer (3 months, reduced area)")\nprint(f"\\n🗺️ CRS: EPSG:32648")\nprint(f" Resolution: 20m (reduced from 10m for smaller data size)")\n\nprint("="*70)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[9]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# 📡 LOAD SENTINEL-2 FROM MICROSOFT PLANETARY COMPUTER\nprint("="*70)\nprint("LOADING SENTINEL-2 L2A")\nprint("="*70)\n\nprint("\\n🔍 Searching for Sentinel-2 scenes...")\nquery_s2 = catalog.search(\n collections=["sentinel-2-l2a"],\n bbox=bbox,\n datetime=time_range,\n query={"eo:cloud_cover": {"lt": 30}} # Cloud cover < 30% (giảm từ 50%)\n)\n\nitems_s2 = list(query_s2.item_collection())\nprint(f"✅ Found {len(items_s2)} Sentinel-2 scenes")\n\n# GIỚI HẠN SỐ LƯỢNG SCENES cho PC cá nhân\nmax_scenes = 12 # Giảm xuống 12 scenes để tối ưu cho PC\nif len(items_s2) > max_scenes:\n print(f"⚠️ Limiting to {max_scenes} scenes for personal computer")\n # Chọn scenes đều đặn trong khoảng thời gian\n step = len(items_s2) // max_scenes\n items_s2 = items_s2[::step][:max_scenes]\n print(f" Selected {len(items_s2)} scenes evenly distributed")\n\nif len(items_s2) > 0:\n # Show first few scenes\n print(f"\\n📋 Sample scenes:")\n for i, item in enumerate(items_s2[:5]):\n date = item.datetime.strftime("%Y-%m-%d")\n cloud = item.properties.get("eo:cloud_cover", "N/A")\n print(f" [{i+1}] {date} - Cloud: {cloud}%")\n \n # Re-sign items to ensure fresh URLs (keep as pystac objects)\n print(f"\\n🔑 Signing STAC items...")\n items_s2 = [planetary_computer.sign(item) for item in items_s2]\n \n # Load Sentinel-2 data (without Dask chunks)\n print(f"\\n⏳ Loading Sentinel-2 data...")\n ds_s2 = stac_load(\n items_s2,\n bands=["B04", "B08", "SCL"], # Red (B04), NIR (B08), Scene Classification (SCL)\n crs="EPSG:32648",\n resolution=20, # 20m resolution (4x smaller data than 10m)\n bbox=bbox,\n patch_url=planetary_computer.sign, # Re-sign URLs during loading\n fail_on_error=False, # Skip problematic tiles instead of crashing\n )\n \n # Rename bands to simpler names\n ds_s2 = ds_s2.rename({"B04": "red", "B08": "nir", "SCL": "scl"})\n \n print(f"\\n✅ Sentinel-2 loaded!")\n print(f" Shape: {dict(ds_s2.dims)}")\n print(f" Variables: {list(ds_s2.data_vars)}")\n display(ds_s2)\nelse:\n print(f"❌ No Sentinel-2 scenes found")\n\n ds_s2 = Noneprint("="*70)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[10]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# 📡 LOAD SENTINEL-1 FROM MICROSOFT PLANETARY COMPUTER\nprint("="*70)\nprint("LOADING SENTINEL-1 RTC")\nprint("="*70)\n\nprint("\\n🔍 Searching for Sentinel-1 scenes...")\nquery_s1 = catalog.search(\n collections=["sentinel-1-rtc"],\n bbox=bbox,\n datetime=time_range,\n)\n\nitems_s1 = list(query_s1.item_collection())\nprint(f"✅ Found {len(items_s1)} Sentinel-1 scenes")\n\n# GIỚI HẠN SỐ LƯỢNG SCENES cho PC cá nhân\nmax_scenes = 12 # Giảm xuống 12 scenes để tối ưu cho PC\nif len(items_s1) > max_scenes:\n print(f"⚠️ Limiting to {max_scenes} scenes for personal computer")\n # Chọn scenes đều đặn trong khoảng thời gian\n step = len(items_s1) // max_scenes\n items_s1 = items_s1[::step][:max_scenes]\n print(f" Selected {len(items_s1)} scenes evenly distributed")\n\nif len(items_s1) > 0:\n # Show first few scenes\n print(f"\\n📋 Sample scenes:")\n for i, item in enumerate(items_s1[:5]):\n date = item.datetime.strftime("%Y-%m-%d")\n orbit = item.properties.get("sat:orbit_state", "N/A")\n print(f" [{i+1}] {date} - Orbit: {orbit}")\n \n # Re-sign items to ensure fresh URLs (keep as pystac objects)\n print(f"\\n🔑 Signing STAC items...")\n items_s1 = [planetary_computer.sign(item) for item in items_s1]\n \n # Load Sentinel-1 data (without Dask chunks)\n print(f"\\n⏳ Loading Sentinel-1 data...")\n ds_s1 = stac_load(\n items_s1,\n bands=["vv", "vh"], # VV and VH polarizations\n crs="EPSG:32648",\n resolution=20, # 20m resolution (4x smaller data than 10m)\n bbox=bbox,\n patch_url=planetary_computer.sign, # Re-sign URLs during loading\n fail_on_error=False, # Skip problematic tiles instead of crashing\n )\n \n # Convert to dB (Microsoft S1 is in linear power)\n print(f"\\n🔄 Converting to dB...")\n ds_s1[\'vv_db\'] = 10 * np.log10(ds_s1[\'vv\'].where(ds_s1[\'vv\'] > 0))\n ds_s1[\'vh_db\'] = 10 * np.log10(ds_s1[\'vh\'].where(ds_s1[\'vh\'] > 0))\n \n print(f"\\n✅ Sentinel-1 loaded!")\n print(f" Shape: {dict(ds_s1.dims)}")\n print(f" Variables: {list(ds_s1.data_vars)}")\n display(ds_s1)\nelse:\n print(f"❌ No Sentinel-1 scenes found")\n\n ds_s1 = Noneprint("="*70)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[11]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# 🌿 CALCULATE NDVI AND PROCESS DATA\nprint("="*70)\nprint("DATA PROCESSING")\nprint("="*70)\n\nif ds_s2 is not None:\n print("\\n[1] Calculating NDVI...")\n # NDVI = (NIR - Red) / (NIR + Red)\n ndvi = (ds_s2[\'nir\'] - ds_s2[\'red\']) / (ds_s2[\'nir\'] + ds_s2[\'red\'] + 1e-8)\n \n print(f"✅ NDVI calculated")\n print(f" Shape: {ndvi.shape}")\n print(f" Time steps: {len(ndvi.time)}")\n \n # Cloud masking using SCL band\n print(f"\\n[2] Applying cloud mask...")\n # SCL values: 1=defective, 3=cloud shadow, 8=cloud medium, 9=cloud high, 10=cirrus\n cloud_mask = ds_s2[\'scl\'].isin([1, 3, 8, 9, 10])\n ndvi_masked = ndvi.where(~cloud_mask)\n \n print(f"✅ Cloud mask applied")\n \n # Temporal aggregation (mean over time)\n print(f"\\n[3] Computing mean NDVI across time...")\n ndvi_mean = ndvi_masked.mean(dim=\'time\')\n \n # Data already in memory, no need to compute() again\n print(f"✅ Mean NDVI computed")\n print(f" Shape: {ndvi_mean.shape}")\n \nelse:\n print("❌ No Sentinel-2 data to process")\n ndvi_mean = None\n\nprint("="*70)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# 🎯 EXTRACT TRAINING DATA FEATURES\nprint("="*70)\nprint("FEATURE EXTRACTION")\nprint("="*70)\n\n# Check if required data is available\nif \'ndvi_mean\' not in globals() or \'ds_s1\' not in globals():\n print("❌ Error: Please run Cell 6 (DATA PROCESSING) first!")\n print(" Required variables: ndvi_mean, ds_s1")\n raise RuntimeError("Missing required data. Run cells in order: Cell 4 → Cell 5 → Cell 6 → Cell 7")\n\n# Load training shapefile\nimport geopandas as gpd\n\ntrain_path = \'train/ST_training data_updated_1130points_new.shp\'\nprint(f"\\n[1] Loading training data from: {train_path}")\ntrain_gdf = gpd.read_file(train_path)\n\n# Ensure CRS matches\nif train_gdf.crs != \'EPSG:32648\':\n print(f" Reprojecting from {train_gdf.crs} to EPSG:32648...")\n train_gdf = train_gdf.to_crs(\'EPSG:32648\')\n\nprint(f"✅ Loaded {len(train_gdf)} training points")\nprint(f" Available columns: {list(train_gdf.columns)}")\n\n# Auto-detect label column (look for common names)\nlabel_column = None\nfor col in [\'HT_code\', \'Ma_LU\', \'LU2022\', \'class\', \'Class\', \'CLASS\', \'label\', \'Label\', \'LABEL\', \'LU_CODE\', \'LU_code\']:\n if col in train_gdf.columns:\n label_column = col\n break\n\nif label_column is None:\n print(f"❌ Cannot find label column. Available columns: {list(train_gdf.columns)}")\n print(f" Please check your shapefile and update the code.")\nelse:\n print(f" Using label column: \'{label_column}\'")\n print(f" Classes: {sorted(train_gdf[label_column].unique())}")\n \n # Extract features at each training point\n print(f"\\n[2] Extracting features at training points...")\n \n features = []\n labels = []\n skipped = 0\n \n for idx, row in train_gdf.iterrows():\n point = row.geometrychro\n x_coord = point.x\n y_coord = point.y\n label = row[label_column]\n \n # Extract NDVI at this location\n if ndvi_mean is not None and ds_s1 is not None:\n try:\n ndvi_val = ndvi_mean.sel(x=x_coord, y=y_coord, method=\'nearest\').values\n \n # Extract Sentinel-1 VH/VV at this location (mean across time)\n # Data already in memory, no need to compute()\n vh_val = ds_s1[\'vh_db\'].sel(x=x_coord, y=y_coord, method=\'nearest\').mean(dim=\'time\').values\n vv_val = ds_s1[\'vv_db\'].sel(x=x_coord, y=y_coord, method=\'nearest\').mean(dim=\'time\').values\n \n # Create feature vector: [NDVI, VH_dB, VV_dB]\n feature_vec = [ndvi_val, vh_val, vv_val]\n \n # Only add if all features are valid (not NaN)\n if not np.isnan(feature_vec).any():\n features.append(feature_vec)\n labels.append(label)\n else:\n skipped += 1\n except Exception as e:\n # Skip points outside the data extent\n skipped += 1\n continue\n \n features = np.array(features)\n labels = np.array(labels)\n \n print(f"✅ Extracted features for {len(features)} valid points")\n print(f" Skipped {skipped} points (outside extent or NaN values)")\n print(f" Feature shape: {features.shape}")\n print(f" Feature names: [\'NDVI_mean\', \'VH_dB_mean\', \'VV_dB_mean\']")\n print(f"\\n Class distribution:")\n unique, counts = np.unique(labels, return_counts=True)\n for cls, cnt in zip(unique, counts):\n print(f" Class {cls}: {cnt} samples ({cnt/len(labels)*100:.1f}%)")\n\nprint("="*70)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[21]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# 🤖 TRAIN XGBOOST MODEL ON GPU (RTX 4060)\nprint("="*70)\nprint("MODEL TRAINING - GPU ACCELERATED")\nprint("="*70)\n\nfrom xgboost import XGBClassifier\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.preprocessing import LabelEncoder\nfrom sklearn.metrics import classification_report, confusion_matrix, ConfusionMatrixDisplay\nimport matplotlib.pyplot as plt\n\n# Encode labels to ensure they are 0, 1, 2, ... n-1\nprint("\\n[1] Encoding labels...")\nlabel_encoder = LabelEncoder()\nlabels_encoded = label_encoder.fit_transform(labels)\nprint(f"✅ Original classes: {label_encoder.classes_}")\nprint(f" Encoded as: {np.unique(labels_encoded)}")\n\n# Split data\nprint("\\n[2] Splitting data (80% train, 20% test)...")\nX_train, X_test, y_train, y_test = train_test_split(\n features, labels_encoded, test_size=0.2, random_state=42, stratify=labels_encoded\n)\nprint(f"✅ Training samples: {len(X_train)}")\nprint(f" Testing samples: {len(X_test)}")\n\n# Train XGBoost on GPU\nprint("\\n[3] Training XGBoost classifier on RTX 4060 GPU...")\nprint(" GPU Settings: device=\'cuda:0\'")\n\nxgb_model = XGBClassifier(\n n_estimators=100,\n max_depth=20,\n learning_rate=0.1,\n device=\'cuda:0\', # Use GPU (updated from deprecated gpu_id)\n tree_method=\'hist\', # Use hist with device for GPU training\n random_state=42,\n eval_metric=\'mlogloss\', # Multi-class log loss\n verbosity=1 # Show GPU training progress\n)\n\nxgb_model.fit(X_train, y_train)\nprint(f"✅ Model trained on GPU")\n\n# Evaluate\nprint("\\n[4] Evaluating model...")\ntrain_score = xgb_model.score(X_train, y_train)\ntest_score = xgb_model.score(X_test, y_test)\nprint(f"✅ Training accuracy: {train_score:.4f}")\nprint(f" Testing accuracy: {test_score:.4f}")\n\n# Classification report\nprint("\\n[5] Classification Report:")\ny_pred = xgb_model.predict(X_test)\nprint(classification_report(y_test, y_pred, target_names=[str(c) for c in label_encoder.classes_]))\n\n# Confusion matrix\nprint("\\n[6] Confusion Matrix:")\nfig, ax = plt.subplots(figsize=(10, 8))\ncm = confusion_matrix(y_test, y_pred)\ndisp = ConfusionMatrixDisplay(confusion_matrix=cm, display_labels=label_encoder.classes_)\ndisp.plot(ax=ax, cmap=\'Blues\', values_format=\'d\')\nplt.title(\'Confusion Matrix - XGBoost GPU Model (RTX 4060)\')\nplt.tight_layout()\nplt.show()\n\nprint("="*70)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[23]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# 💾 SAVE MODEL AND CLEANUP\nprint("="*70)\nprint("SAVING MODEL & CLEANUP")\nprint("="*70)\n\nimport joblib\nfrom datetime import datetime\n\n# Save model and label encoder\nmodel_filename = f"model_train/model_xgboost_gpu_{datetime.now().strftime(\'%Y%m%d_%H%M%S\')}.joblib"\nprint(f"\\n[1] Saving model to: {model_filename}")\njoblib.dump({\'model\': xgb_model, \'label_encoder\': label_encoder}, model_filename)\nprint(f"✅ Model and label encoder saved")\n\n# Save model info\ninfo = {\n "timestamp": datetime.now().isoformat(),\n "data_source": "Microsoft Planetary Computer STAC",\n "collections": ["sentinel-2-l2a", "sentinel-1-rtc"],\n "features": ["NDVI_mean", "VH_dB_mean", "VV_dB_mean"],\n "training_samples": len(X_train),\n "testing_samples": len(X_test),\n "train_accuracy": float(train_score),\n "test_accuracy": float(test_score),\n "model_type": "XGBClassifier",\n "device": "cuda:0",\n "gpu_device": "RTX 4060",\n "tree_method": "hist",\n "n_estimators": 100,\n "max_depth": 20,\n "learning_rate": 0.1\n}\n\nimport json\ninfo_filename = model_filename.replace(\'.joblib\', \'_info.json\')\nwith open(info_filename, \'w\') as f:\n json.dump(info, f, indent=2)\nprint(f"✅ Model info saved to: {info_filename}")\n\n# No cleanup needed (Dask removed)\nprint("\\n[2] Cleanup complete")\n\nprint("="*70)\n\nprint("\\n" + "="*70)\n\nprint("🎉 TRAINING COMPLETE!")\n')
|
|
||||||
|
|
||||||
File diff suppressed because one or more lines are too long
@@ -1,191 +0,0 @@
|
|||||||
#!/usr/bin/env python
|
|
||||||
# coding: utf-8
|
|
||||||
|
|
||||||
# In[1]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '%matplotlib inline\n\nimport importlib\nimport new_import_ODC \n\nimportlib.reload(new_import_ODC)\n\nfrom new_import_ODC import *\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[2]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Dask gateway\ncluster, client = notebook_utils.initialize_dask(use_gateway=True, workers=(1,4))\ndc = datacube.Datacube()\n\n# Configure s3 access\nconfigure_s3_access(aws_unsigned=False, requester_pays=True, client=client)\n\nclient\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[3]:
|
|
||||||
|
|
||||||
|
|
||||||
## cấu hình thời gian lấy ảnh và tọa độ
|
|
||||||
date_range = ('2022-09-01', '2023-10-01')
|
|
||||||
longtitude_range = (105.86575, 105.94120)
|
|
||||||
latitude_range = (9.65070, 9.69850)
|
|
||||||
|
|
||||||
|
|
||||||
# In[4]:
|
|
||||||
|
|
||||||
|
|
||||||
## truy vấn ảnh vệ tinh sen2
|
|
||||||
data = load_data(dc, date_range, longtitude_range, latitude_range)
|
|
||||||
notebook_utils.heading(notebook_utils.xarray_object_size(data))
|
|
||||||
display(data)
|
|
||||||
|
|
||||||
|
|
||||||
# In[5]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Tiến hành loại bỏ các vị trí bị mây ảnh hưởng\nresult = mask_clean(data)\nprogress(result)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[6]:
|
|
||||||
|
|
||||||
|
|
||||||
# Tiến hành tính toán NDVI
|
|
||||||
ds1 = calculate_indices(result, index='NDVI', satellite_mission='s2')
|
|
||||||
ndvi = ds1["NDVI"]
|
|
||||||
display(ndvi)
|
|
||||||
|
|
||||||
|
|
||||||
# In[7]:
|
|
||||||
|
|
||||||
|
|
||||||
## ảnh NDVI chưa điền mây (fill nan)
|
|
||||||
plt.imshow(ndvi.isel(time=50))
|
|
||||||
|
|
||||||
|
|
||||||
# In[8]:
|
|
||||||
|
|
||||||
|
|
||||||
# đặt thời gian các mùa
|
|
||||||
time_split = [slice('2022-09-01', '2023-01-01'),
|
|
||||||
slice('2023-01-01', '2023-05-01'),
|
|
||||||
slice('2023-05-01', '2023-07-01'),
|
|
||||||
slice('2023-07-01', '2023-10-01')]
|
|
||||||
|
|
||||||
# Điền mây ở các vị trí mang giá trị nan (fill nan)
|
|
||||||
fill_nan_ndvi = fill_nan(ndvi, time_split)
|
|
||||||
|
|
||||||
# In kết quả ảnh ndvi đã điền mây (đã fill nan)
|
|
||||||
plt.imshow(fill_nan_ndvi.isel(time=50))
|
|
||||||
|
|
||||||
|
|
||||||
# In[9]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', "## tính ndvi theo tháng\naverage_ndvi = fill_nan_ndvi.resample(time='1M').mean().persist()\nprogress(average_ndvi)\n")
|
|
||||||
|
|
||||||
|
|
||||||
# In[10]:
|
|
||||||
|
|
||||||
|
|
||||||
# compute average_ndvi
|
|
||||||
average_ndvi = average_ndvi.compute()
|
|
||||||
|
|
||||||
|
|
||||||
# In[11]:
|
|
||||||
|
|
||||||
|
|
||||||
# load dữ liệu sen1
|
|
||||||
coordinates = (longtitude_range, latitude_range)
|
|
||||||
dsvh, dsvv = load_data_sen1(dc, date_range, coordinates)
|
|
||||||
average_vv = calculate_average(dsvv, time_pattern='1M')
|
|
||||||
average_vh = calculate_average(dsvh, time_pattern='1M')
|
|
||||||
|
|
||||||
|
|
||||||
# In[12]:
|
|
||||||
|
|
||||||
|
|
||||||
# load model RF
|
|
||||||
loaded_model = joblib.load(os.path.join("model_train", "model_odc.joblib"))
|
|
||||||
|
|
||||||
# dự đoán
|
|
||||||
data_array = predict(loaded_model, data.rio.crs, average_ndvi, average_vh, average_vv)
|
|
||||||
|
|
||||||
|
|
||||||
# In[13]:
|
|
||||||
|
|
||||||
|
|
||||||
# cấu hình màu cho các loại đất
|
|
||||||
colors = [
|
|
||||||
"#abcee9",
|
|
||||||
"#ffef44",
|
|
||||||
"#c4ff9e",
|
|
||||||
"#ffd6a8",
|
|
||||||
"#93ddda",
|
|
||||||
"#1aeef7",
|
|
||||||
"#ffa7f2",
|
|
||||||
"#33ee33"
|
|
||||||
]
|
|
||||||
labels = [
|
|
||||||
"Lúa tôm",
|
|
||||||
"Lúa",
|
|
||||||
"CHN",
|
|
||||||
"CLN",
|
|
||||||
"TS",
|
|
||||||
"Sông",
|
|
||||||
"Đất xây dựng",
|
|
||||||
"Rừng"
|
|
||||||
]
|
|
||||||
# hiển thị phân loại sử dụng đất
|
|
||||||
cmap = ListedColormap(colors)
|
|
||||||
img = data_array.plot(cmap=cmap, add_colorbar=False)
|
|
||||||
cbar = plt.colorbar(img)
|
|
||||||
cbar.ax.set_yticklabels(labels)
|
|
||||||
plt.title("Phân loại sử dụng đất")
|
|
||||||
plt.axis('off')
|
|
||||||
plt.show()
|
|
||||||
|
|
||||||
|
|
||||||
# In[14]:
|
|
||||||
|
|
||||||
|
|
||||||
## cấu hình shapefile ranh giới thuận hòa và vh vv file
|
|
||||||
thuanhoa_path = "ThuanHoa/region/ST_ThuanHoa_Boundaryofficially.shp"
|
|
||||||
|
|
||||||
# cắt theo ranh giới xã thuận hòa
|
|
||||||
region_result = cut_according_shp(thuanhoa_path, average_ndvi, data_array)
|
|
||||||
|
|
||||||
|
|
||||||
# In[15]:
|
|
||||||
|
|
||||||
|
|
||||||
# hiển thị kết quả phân loại sử dụng đất
|
|
||||||
colorval = list(range(len(colors)))
|
|
||||||
options = {
|
|
||||||
'title': 'Phân loại sử dụng đất',
|
|
||||||
'cmap': colors,
|
|
||||||
'clim': (0, 8),
|
|
||||||
'aspect': 'equal',
|
|
||||||
'colorbar_opts': {
|
|
||||||
'major_label_overrides': dict(zip(colorval, labels)),
|
|
||||||
'major_label_text_align': 'left',
|
|
||||||
'ticker': FixedTicker(ticks=colorval),
|
|
||||||
},
|
|
||||||
}
|
|
||||||
|
|
||||||
region_result.hvplot(
|
|
||||||
rasterize = True, # Use Datashader, particularly useful for dask arrays
|
|
||||||
aggregator = reductions.mode(), # Datashader selects mode value, requires 'hv.Image'
|
|
||||||
).options(opts.Image(**options))
|
|
||||||
|
|
||||||
|
|
||||||
# In[16]:
|
|
||||||
|
|
||||||
|
|
||||||
# Lưu lại kết quả
|
|
||||||
region_result.rio.to_raster("KetQuaPhanLoaiDatODC.tif")
|
|
||||||
|
|
||||||
|
|
||||||
# In[17]:
|
|
||||||
|
|
||||||
|
|
||||||
# đóng client, cluster
|
|
||||||
client.close()
|
|
||||||
cluster.close()
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
@@ -1,117 +0,0 @@
|
|||||||
#!/usr/bin/env python
|
|
||||||
# coding: utf-8
|
|
||||||
|
|
||||||
# In[1]:
|
|
||||||
|
|
||||||
|
|
||||||
# Khai báo các thư viện cần thiết
|
|
||||||
from new_import_ODC import *
|
|
||||||
|
|
||||||
# Khai báo đường dẫn đến kết quả phân loại và dữ liệu của địa phương
|
|
||||||
KD_path = "ThuanHoa/KhoanhDat/ThuanHoa_TKDD2022.shp"
|
|
||||||
KetQuaPhanLoaiDat = "KetQuaPhanLoaiDatODC.tif"
|
|
||||||
|
|
||||||
|
|
||||||
# In[2]:
|
|
||||||
|
|
||||||
|
|
||||||
# khai báo các loại đất từ dữ liệu kiểm kê ứng với các hiện trạng được phân loại từ viễn thám
|
|
||||||
CODE_MAP = {
|
|
||||||
"BHK": 2,
|
|
||||||
"CLN": 3,
|
|
||||||
"DGD": 6,
|
|
||||||
"DGT": 6,
|
|
||||||
"DNL": 6,
|
|
||||||
"DRA": 6,
|
|
||||||
"DSH": 6,
|
|
||||||
"DTL": 5,
|
|
||||||
"DTS": 6,
|
|
||||||
"DYT": 6,
|
|
||||||
"LUC": 1,
|
|
||||||
"NKH": 3,
|
|
||||||
"NTD": 6,
|
|
||||||
"NTS": 4,
|
|
||||||
"ONT": 6,
|
|
||||||
"SKC": 6,
|
|
||||||
"SKX": 6,
|
|
||||||
"SON": 5,
|
|
||||||
"TMD": 6,
|
|
||||||
"TON": 6,
|
|
||||||
"TSC": 6,
|
|
||||||
}
|
|
||||||
|
|
||||||
# Khai báo các nhãn phân loại đất ứng với 3 loại đất chính
|
|
||||||
HT_MAP = {
|
|
||||||
"NN": {"name": "Đất Nông Nghiệp", "data": [1, 2, 3, 4]},
|
|
||||||
"PNN": {"name": "Đất Phi Nông Nghiệp", "data": [6]},
|
|
||||||
"TQ": {"name": "Đất Thổ Quả", "data": [15]},
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
# In[3]:
|
|
||||||
|
|
||||||
|
|
||||||
# Tiến hành chồng lắp
|
|
||||||
result = compare(KD_path, KetQuaPhanLoaiDat, CODE_MAP, HT_MAP)
|
|
||||||
|
|
||||||
|
|
||||||
# In[4]:
|
|
||||||
|
|
||||||
|
|
||||||
# cấu hình màu cho các loại sử dụng đất
|
|
||||||
colors = [
|
|
||||||
"#abcee9",
|
|
||||||
"#ffffc0",
|
|
||||||
"#c4ff9e",
|
|
||||||
"#ffd6a8",
|
|
||||||
"#93ddda",
|
|
||||||
"#1aeef7",
|
|
||||||
"#ffa7f2",
|
|
||||||
"#33ee33",
|
|
||||||
]
|
|
||||||
labels = ["Lúa tôm", "Lúa", "CHN", "CLN", "TS", "Sông", "Đất xây dựng", "Rừng"]
|
|
||||||
|
|
||||||
|
|
||||||
# In[5]:
|
|
||||||
|
|
||||||
|
|
||||||
# Lưu kết quả
|
|
||||||
save_result(result, HT_MAP)
|
|
||||||
|
|
||||||
|
|
||||||
# In[6]:
|
|
||||||
|
|
||||||
|
|
||||||
# hiển thị kết quả
|
|
||||||
xx = []
|
|
||||||
|
|
||||||
for k, v in result.items():
|
|
||||||
rs = merge_arrays(v, nodata=np.nan)
|
|
||||||
xx.append(rs.squeeze(drop=True))
|
|
||||||
xx = xr.concat(xx, pd.Index([HT_MAP[x]["name"] for x in HT_MAP], name="name"))
|
|
||||||
|
|
||||||
colorval = list(range(len(colors)))
|
|
||||||
options = {
|
|
||||||
"cmap": colors,
|
|
||||||
"clim": (0, 8),
|
|
||||||
"aspect": "equal",
|
|
||||||
"height": 400,
|
|
||||||
"colorbar_opts": {
|
|
||||||
"major_label_overrides": dict(zip(colorval, labels)),
|
|
||||||
"major_label_text_align": "left",
|
|
||||||
"ticker": FixedTicker(ticks=colorval),
|
|
||||||
},
|
|
||||||
}
|
|
||||||
|
|
||||||
xx.hvplot(
|
|
||||||
groupby="name",
|
|
||||||
rasterize=True, # Use Datashader, particularly useful for dask arrays
|
|
||||||
aggregator=reductions.mode(), # Datashader selects mode value, requires 'hv.Image'
|
|
||||||
).options(opts.Image(**options))
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
File diff suppressed because one or more lines are too long
@@ -1,18 +0,0 @@
|
|||||||
import joblib
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
cache_file = "dataset_cache/training_data_2d.joblib"
|
|
||||||
data = joblib.load(cache_file)
|
|
||||||
X = np.array(data['X'])
|
|
||||||
y = np.array(data['y'])
|
|
||||||
|
|
||||||
print("X shape:", X.shape)
|
|
||||||
print("X mean:", np.mean(X))
|
|
||||||
print("X std:", np.std(X))
|
|
||||||
print("X min:", np.min(X))
|
|
||||||
print("X max:", np.max(X))
|
|
||||||
print("Any NaN:", np.isnan(X).any())
|
|
||||||
|
|
||||||
for i in range(6):
|
|
||||||
print(f"Channel {i} mean: {np.mean(X[:, i, :, :]):.4f}, min: {np.min(X[:, i, :, :]):.4f}, max: {np.max(X[:, i, :, :]):.4f}")
|
|
||||||
|
|
||||||
@@ -1,10 +0,0 @@
|
|||||||
import joblib
|
|
||||||
import geopandas as gpd
|
|
||||||
from shapely.geometry import Point
|
|
||||||
|
|
||||||
data = joblib.load('dataset_cache/training_data_2d.joblib')
|
|
||||||
X, y = data['X'], data['y']
|
|
||||||
print(f"X shape: {X.shape}, y shape: {y.shape}")
|
|
||||||
|
|
||||||
gdf = gpd.read_file("train/ST_training_data_updated_1130points_new.shp")
|
|
||||||
print("Total points:", len(gdf))
|
|
||||||
@@ -1,7 +0,0 @@
|
|||||||
import joblib
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
data = joblib.load('dataset_cache/training_data.joblib')
|
|
||||||
X, y = data['X'], data['y']
|
|
||||||
print(f"X shape: {X.shape}")
|
|
||||||
print(f"y shape: {y.shape}")
|
|
||||||
@@ -1,8 +0,0 @@
|
|||||||
import joblib
|
|
||||||
import numpy as np
|
|
||||||
cache_file = "dataset_cache/training_data_2d.joblib"
|
|
||||||
data = joblib.load(cache_file)
|
|
||||||
X = np.array(data['X'])
|
|
||||||
b2 = X[:, 0, :, :]
|
|
||||||
print("Zeros in B2:", np.sum(b2 == 0) / b2.size)
|
|
||||||
print("X shape:", X.shape)
|
|
||||||
@@ -1,25 +0,0 @@
|
|||||||
import glob, json
|
|
||||||
|
|
||||||
changed_files = []
|
|
||||||
for file_path in glob.glob('*.ipynb'):
|
|
||||||
with open(file_path, 'r', encoding='utf-8') as f:
|
|
||||||
nb = json.load(f)
|
|
||||||
|
|
||||||
changed = False
|
|
||||||
for cell in nb.get('cells', []):
|
|
||||||
if cell.get('cell_type') == 'code':
|
|
||||||
source = cell.get('source', [])
|
|
||||||
for i, line in enumerate(source):
|
|
||||||
if 'time=50' in line:
|
|
||||||
source[i] = line.replace('time=50', 'time=0')
|
|
||||||
changed = True
|
|
||||||
if 'load_data_sen1(dc,' in line:
|
|
||||||
source[i] = line.replace('load_data_sen1(dc,', 'load_data_sen1(None,')
|
|
||||||
changed = True
|
|
||||||
|
|
||||||
if changed:
|
|
||||||
with open(file_path, 'w', encoding='utf-8') as f:
|
|
||||||
json.dump(nb, f, indent=1)
|
|
||||||
changed_files.append(file_path)
|
|
||||||
|
|
||||||
print('Fixed issues in:', changed_files)
|
|
||||||
@@ -1,19 +0,0 @@
|
|||||||
import joblib
|
|
||||||
import geopandas as gpd
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
cache_file = "dataset_cache/training_data_2d.joblib"
|
|
||||||
data = joblib.load(cache_file)
|
|
||||||
X = data['X']
|
|
||||||
print("X shape:", len(X))
|
|
||||||
|
|
||||||
gdf = gpd.read_file("train/ST_training_data_updated_1130points_new.shp")
|
|
||||||
gdf = gdf.to_crs("EPSG:32648")
|
|
||||||
print("gdf length:", len(gdf))
|
|
||||||
|
|
||||||
if len(X) == len(gdf):
|
|
||||||
y = [(row['HT_code'] - 1) for idx, row in gdf.iterrows()]
|
|
||||||
joblib.dump({'X': X, 'y': y}, cache_file)
|
|
||||||
print("Fixed y in cache! Saved.")
|
|
||||||
else:
|
|
||||||
print("Lengths do not match, cannot fix automatically.")
|
|
||||||
@@ -1,20 +0,0 @@
|
|||||||
import json
|
|
||||||
|
|
||||||
def fix_import(file_path):
|
|
||||||
with open(file_path, 'r', encoding='utf-8') as f:
|
|
||||||
nb = json.load(f)
|
|
||||||
changed = False
|
|
||||||
for cell in nb.get('cells', []):
|
|
||||||
if cell.get('cell_type') == 'code':
|
|
||||||
source = cell.get('source', [])
|
|
||||||
if isinstance(source, list):
|
|
||||||
for i, line in enumerate(source):
|
|
||||||
if "from new_import import *" in line:
|
|
||||||
source[i] = line.replace("from new_import import *", "from new_import_ODC import *")
|
|
||||||
changed = True
|
|
||||||
if changed:
|
|
||||||
with open(file_path, 'w', encoding='utf-8') as f:
|
|
||||||
json.dump(nb, f, indent=1)
|
|
||||||
print(f"Fixed {file_path}")
|
|
||||||
|
|
||||||
fix_import('new_train.ipynb')
|
|
||||||
@@ -1,20 +0,0 @@
|
|||||||
import re
|
|
||||||
|
|
||||||
filepath = "01.train_ODC.py"
|
|
||||||
with open(filepath, 'r', encoding='utf-8') as f:
|
|
||||||
content = f.read()
|
|
||||||
|
|
||||||
# Find the run_cell_magic line
|
|
||||||
pattern = re.compile(r"get_ipython\(\)\.run_cell_magic\('time', '', '(# 🤖 LAND USE CLASSIFICATION MODEL TRAINING.*?)(?=\n')\n'", re.DOTALL)
|
|
||||||
|
|
||||||
def repl(match):
|
|
||||||
# Get the inner string and escape all actual newlines with \n
|
|
||||||
inner = match.group(1)
|
|
||||||
inner = inner.replace('\n', '\\n')
|
|
||||||
return f"get_ipython().run_cell_magic('time', '', '{inner}')"
|
|
||||||
|
|
||||||
content = pattern.sub(repl, content)
|
|
||||||
|
|
||||||
with open(filepath, 'w', encoding='utf-8') as f:
|
|
||||||
f.write(content)
|
|
||||||
print("Fixed 01.train_ODC.py syntax")
|
|
||||||
@@ -1,22 +0,0 @@
|
|||||||
import json
|
|
||||||
|
|
||||||
def fix_filename(file_path):
|
|
||||||
with open(file_path, 'r', encoding='utf-8') as f:
|
|
||||||
nb = json.load(f)
|
|
||||||
changed = False
|
|
||||||
for cell in nb.get('cells', []):
|
|
||||||
if cell.get('cell_type') == 'code':
|
|
||||||
source = cell.get('source', [])
|
|
||||||
if isinstance(source, list):
|
|
||||||
for i, line in enumerate(source):
|
|
||||||
if "ST_training data_updated_1130points.shp" in line:
|
|
||||||
source[i] = line.replace("ST_training data_updated_1130points.shp", "ST_training_data_updated_1130points.shp")
|
|
||||||
changed = True
|
|
||||||
if changed:
|
|
||||||
with open(file_path, 'w', encoding='utf-8') as f:
|
|
||||||
json.dump(nb, f, indent=1)
|
|
||||||
print(f"Fixed typo in {file_path}")
|
|
||||||
|
|
||||||
import glob
|
|
||||||
for nb in glob.glob("*.ipynb"):
|
|
||||||
fix_filename(nb)
|
|
||||||
@@ -1,58 +0,0 @@
|
|||||||
import nbformat as nbf
|
|
||||||
|
|
||||||
nb = nbf.v4.new_notebook()
|
|
||||||
|
|
||||||
text_1 = """# Script Tải Dữ liệu Vệ tinh (Cache) qua Google Colab
|
|
||||||
Mục đích của Notebook này là mượn sức mạnh đường truyền và RAM của Google Colab để tải 270 ảnh Sentinel-2 & Sentinel-1 từ Microsoft Planetary Computer. Sau khi xử lý nội suy, nó sẽ sinh ra một file cache `.joblib` duy nhất chứa toàn bộ mảng dữ liệu.
|
|
||||||
Bạn chỉ cần tải file `.joblib` đó về máy là xong!"""
|
|
||||||
|
|
||||||
code_1 = """!pip install planetary-computer pystac-client odc-stac geopandas rasterio xarray joblib scikit-learn xgboost lightgbm"""
|
|
||||||
|
|
||||||
code_2 = """from google.colab import drive
|
|
||||||
drive.mount('/content/drive')"""
|
|
||||||
|
|
||||||
text_2 = """## Hướng dẫn:
|
|
||||||
1. Nén toàn bộ thư mục `remote-sensing` ở máy tính của bạn thành file `remote-sensing.zip`.
|
|
||||||
2. Upload file `remote-sensing.zip` đó lên Google Drive (để ngay ngoài cùng).
|
|
||||||
3. Chạy ô lệnh bên dưới để giải nén và chuyển vào thư mục dự án."""
|
|
||||||
|
|
||||||
code_3 = """import os
|
|
||||||
import shutil
|
|
||||||
|
|
||||||
# Giải nén dự án từ Google Drive
|
|
||||||
!unzip -q /content/drive/MyDrive/remote-sensing.zip -d /content/
|
|
||||||
os.chdir('/content/remote-sensing')
|
|
||||||
!ls -la"""
|
|
||||||
|
|
||||||
text_3 = """## Bắt đầu tải và Cache Dữ Liệu
|
|
||||||
Chạy một mô hình CPU đơn giản (Decision Tree) để ép hệ thống gọi hàm `FeatureExtractor`. Hàm này sẽ làm mọi việc nặng nhọc: tìm ảnh, ghép mây, tính trung vị và lưu kết quả vào thư mục `dataset_cache/`."""
|
|
||||||
|
|
||||||
code_4 = """# Lệnh này sẽ mất khoảng 5-15 phút để tải toàn bộ ảnh từ Microsoft
|
|
||||||
!python train_land_decision_tree_gpu.py"""
|
|
||||||
|
|
||||||
text_4 = """## Hoàn tất
|
|
||||||
Bạn hãy kiểm tra xem file `.joblib` lớn (khoảng 40-60MB) đã xuất hiện chưa. Nếu rồi, hãy lưu ngược nó lại Google Drive để tải về máy!"""
|
|
||||||
|
|
||||||
code_5 = """# Xem file cache đã được tạo thành công chưa
|
|
||||||
!ls -lh dataset_cache/
|
|
||||||
|
|
||||||
# Copy toàn bộ thư mục cache sang Google Drive để tải về máy dễ dàng
|
|
||||||
!cp -r dataset_cache/ /content/drive/MyDrive/dataset_cache_finished/
|
|
||||||
print("Hoàn thành! Bạn hãy mở Google Drive của mình, tìm thư mục 'dataset_cache_finished' và tải file .joblib mới nhất về máy tính.")"""
|
|
||||||
|
|
||||||
nb['cells'] = [
|
|
||||||
nbf.v4.new_markdown_cell(text_1),
|
|
||||||
nbf.v4.new_code_cell(code_1),
|
|
||||||
nbf.v4.new_code_cell(code_2),
|
|
||||||
nbf.v4.new_markdown_cell(text_2),
|
|
||||||
nbf.v4.new_code_cell(code_3),
|
|
||||||
nbf.v4.new_markdown_cell(text_3),
|
|
||||||
nbf.v4.new_code_cell(code_4),
|
|
||||||
nbf.v4.new_markdown_cell(text_4),
|
|
||||||
nbf.v4.new_code_cell(code_5)
|
|
||||||
]
|
|
||||||
|
|
||||||
with open('Download_Cache_Colab.ipynb', 'w') as f:
|
|
||||||
nbf.write(nb, f)
|
|
||||||
|
|
||||||
print("Created Download_Cache_Colab.ipynb")
|
|
||||||
@@ -1,97 +0,0 @@
|
|||||||
import os
|
|
||||||
import glob
|
|
||||||
import json
|
|
||||||
from tabulate import tabulate
|
|
||||||
|
|
||||||
print("📊 BẢNG SO SÁNH KẾT QUẢ CÁC MÔ HÌNH\n")
|
|
||||||
|
|
||||||
# 1. Phân loại đất
|
|
||||||
print("### 1. Nhóm Phân loại Lớp phủ (Land Classification)")
|
|
||||||
land_data = []
|
|
||||||
|
|
||||||
if os.path.exists("model_xgboost_info.json"):
|
|
||||||
with open("model_xgboost_info.json", 'r') as f:
|
|
||||||
data = json.load(f)
|
|
||||||
params = data.get('params', {})
|
|
||||||
param_str = f"estimators:{params.get('n_estimators')}, depth:{params.get('max_depth')}" if params else "N/A"
|
|
||||||
land_data.append([
|
|
||||||
data.get('model_type', 'XGBoost'),
|
|
||||||
data.get('accuracy', ''),
|
|
||||||
data.get('precision', ''),
|
|
||||||
data.get('recall', ''),
|
|
||||||
data.get('f1_score', ''),
|
|
||||||
param_str
|
|
||||||
])
|
|
||||||
|
|
||||||
for info_file in glob.glob("model_train/*_info.json"):
|
|
||||||
with open(info_file, 'r') as f:
|
|
||||||
data = json.load(f)
|
|
||||||
# Support both 'accuracy' and 'test_accuracy'
|
|
||||||
acc = data.get('accuracy', data.get('test_accuracy', ''))
|
|
||||||
|
|
||||||
f1 = data.get('f1_score', '')
|
|
||||||
precision = data.get('precision', '')
|
|
||||||
recall = data.get('recall', '')
|
|
||||||
|
|
||||||
clf_rep = data.get('classification_report')
|
|
||||||
if isinstance(clf_rep, dict) and 'macro avg' in clf_rep:
|
|
||||||
if not f1:
|
|
||||||
f1 = clf_rep['macro avg'].get('f1-score', '')
|
|
||||||
if not precision:
|
|
||||||
precision = clf_rep['macro avg'].get('precision', '')
|
|
||||||
if not recall:
|
|
||||||
recall = clf_rep['macro avg'].get('recall', '')
|
|
||||||
|
|
||||||
if not acc and not f1:
|
|
||||||
continue
|
|
||||||
|
|
||||||
params = data.get('params', {})
|
|
||||||
param_str = f"estimators:{params.get('n_estimators')}, depth:{params.get('max_depth')}" if params else "N/A"
|
|
||||||
|
|
||||||
if data.get('model_type') == 'RandomForest_RealData':
|
|
||||||
param_str = "estimators:100, depth:15"
|
|
||||||
|
|
||||||
land_data.append([
|
|
||||||
data.get('model_type', ''),
|
|
||||||
acc,
|
|
||||||
precision,
|
|
||||||
recall,
|
|
||||||
f1,
|
|
||||||
param_str
|
|
||||||
])
|
|
||||||
|
|
||||||
if land_data:
|
|
||||||
print(tabulate(land_data, headers=["Model", "Accuracy", "Precision", "Recall", "F1-Score", "Parameters"], tablefmt="github"))
|
|
||||||
print("\n")
|
|
||||||
|
|
||||||
# 2. Xóa mây
|
|
||||||
print("### 2. Nhóm Xóa mây (Cloud Removal)")
|
|
||||||
cloud_data = []
|
|
||||||
for info_file in glob.glob("cloud_removal_model/*_info.json"):
|
|
||||||
with open(info_file, 'r') as f:
|
|
||||||
data = json.load(f)
|
|
||||||
cloud_data.append([
|
|
||||||
data.get('model_type', ''),
|
|
||||||
data.get('epoch', ''),
|
|
||||||
data.get('train_loss', ''),
|
|
||||||
data.get('val_loss', '')
|
|
||||||
])
|
|
||||||
if cloud_data:
|
|
||||||
print(tabulate(cloud_data, headers=["Model", "Epochs", "Train Loss", "Val Loss"], tablefmt="github"))
|
|
||||||
print("\n")
|
|
||||||
|
|
||||||
# 3. Dự báo NDVI
|
|
||||||
print("### 3. Nhóm Dự báo Thực vật (NDVI Forecasting)")
|
|
||||||
ndvi_data = []
|
|
||||||
for info_file in glob.glob("ndvi_forecast_model/*_info.json"):
|
|
||||||
with open(info_file, 'r') as f:
|
|
||||||
data = json.load(f)
|
|
||||||
ndvi_data.append([
|
|
||||||
data.get('model_type', ''),
|
|
||||||
data.get('rmse', ''),
|
|
||||||
data.get('mae', ''),
|
|
||||||
data.get('epoch', 'N/A')
|
|
||||||
])
|
|
||||||
if ndvi_data:
|
|
||||||
print(tabulate(ndvi_data, headers=["Model", "RMSE", "MAE", "Epochs"], tablefmt="github"))
|
|
||||||
print("\n")
|
|
||||||
@@ -1,7 +0,0 @@
|
|||||||
import json
|
|
||||||
nb = json.load(open('01.train_ODC.ipynb'))
|
|
||||||
for idx, cell in enumerate(nb['cells']):
|
|
||||||
if cell['cell_type'] == 'code':
|
|
||||||
print(f"Cell {idx}:")
|
|
||||||
print("".join(cell['source'][:3]))
|
|
||||||
print("-" * 20)
|
|
||||||
@@ -1,18 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU HUẤN LUYỆN MÔ HÌNH CNN (GPU & CACHE)
|
|
||||||
Initializing FeatureExtractor (mode=extended)...
|
|
||||||
📦 Đang load cache: training_data_507bd2ba4ec0d3fe107839cbf73a7a7d.joblib...
|
|
||||||
✅ Loaded 632 samples từ cache!
|
|
||||||
⚡ Đã bỏ qua download위성 data (tiết kiệm thời gian)
|
|
||||||
[CACHE HIT] Using cached dataset with 632 samples
|
|
||||||
Training CNN model...
|
|
||||||
Building CNN model on cuda...
|
|
||||||
Training CNN model with PyTorch...
|
|
||||||
CNN Epoch 10/15, Loss: 1.2171
|
|
||||||
Evaluating model...
|
|
||||||
Generating classification report...
|
|
||||||
Saving model...
|
|
||||||
[MODEL MANAGER] Saving model to: model_train/model_cnn_auto.joblib
|
|
||||||
[MODEL MANAGER] Saving metadata to: model_train/model_cnn_auto_info.json
|
|
||||||
[MODEL MANAGER] Model saved successfully!
|
|
||||||
Training complete!
|
|
||||||
✅ Hoàn thành! Accuracy: 0.5748
|
|
||||||
@@ -1,15 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU HUẤN LUYỆN MÔ HÌNH DECISION TREE (GPU & CACHE)
|
|
||||||
Initializing FeatureExtractor (mode=extended)...
|
|
||||||
📦 Đang load cache: training_data_507bd2ba4ec0d3fe107839cbf73a7a7d.joblib...
|
|
||||||
✅ Loaded 632 samples từ cache!
|
|
||||||
⚡ Đã bỏ qua download위성 data (tiết kiệm thời gian)
|
|
||||||
[CACHE HIT] Using cached dataset with 632 samples
|
|
||||||
Training DECISION_TREE model...
|
|
||||||
Evaluating model...
|
|
||||||
Generating classification report...
|
|
||||||
Saving model...
|
|
||||||
[MODEL MANAGER] Saving model to: model_train/model_decision_tree_auto.joblib
|
|
||||||
[MODEL MANAGER] Saving metadata to: model_train/model_decision_tree_auto_info.json
|
|
||||||
[MODEL MANAGER] Model saved successfully!
|
|
||||||
Training complete!
|
|
||||||
✅ Hoàn thành! Accuracy: 0.5906
|
|
||||||
File diff suppressed because it is too large
Load Diff
File diff suppressed because it is too large
Load Diff
File diff suppressed because it is too large
Load Diff
@@ -1,30 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU HUẤN LUYỆN MÔ HÌNH MOBILENET-LRASPP (GPU & CACHE)
|
|
||||||
Initializing FeatureExtractor (mode=extended)...
|
|
||||||
📦 Đang load cache: training_data_507bd2ba4ec0d3fe107839cbf73a7a7d.joblib...
|
|
||||||
✅ Loaded 632 samples từ cache!
|
|
||||||
⚡ Đã bỏ qua download위성 data (tiết kiệm thời gian)
|
|
||||||
[CACHE HIT] Using cached dataset with 632 samples
|
|
||||||
Training MOBILENET-LRASPP model...
|
|
||||||
Building MobileNetV3 + LR-ASPP model on cuda...
|
|
||||||
[MOBILENET] Class distribution: [ 48 89 3 86 74 38 117 50]
|
|
||||||
[MOBILENET] Class weights: [0.37418982 0.20181024 5.98703525 0.20885013 0.24271772 0.47266082
|
|
||||||
0.15351377 0.35922223]
|
|
||||||
Training MobileNetV3 + LR-ASPP model with PyTorch...
|
|
||||||
MobileNet Epoch 5/25, Train Loss: 1.0614, Val Loss: 1.3584, Val Acc: 40.16%, LR: 0.000800
|
|
||||||
[MOBILENET] Epoch 5/25 - Train Loss: 1.0614, Val Loss: 1.3584, Val Acc: 40.16%
|
|
||||||
MobileNet Epoch 10/25, Train Loss: 0.9054, Val Loss: 0.8582, Val Acc: 59.06%, LR: 0.000800
|
|
||||||
[MOBILENET] Epoch 10/25 - Train Loss: 0.9054, Val Loss: 0.8582, Val Acc: 59.06%
|
|
||||||
MobileNet Epoch 15/25, Train Loss: 0.7728, Val Loss: 0.8231, Val Acc: 61.42%, LR: 0.000800
|
|
||||||
[MOBILENET] Epoch 15/25 - Train Loss: 0.7728, Val Loss: 0.8231, Val Acc: 61.42%
|
|
||||||
MobileNet Epoch 20/25, Train Loss: 0.7125, Val Loss: 0.8783, Val Acc: 59.84%, LR: 0.000400
|
|
||||||
[MOBILENET] Epoch 20/25 - Train Loss: 0.7125, Val Loss: 0.8783, Val Acc: 59.84%
|
|
||||||
[MOBILENET] Early stopping at epoch 24 (best val loss: 0.8029)
|
|
||||||
MobileNet early stopped at epoch 24
|
|
||||||
Evaluating model...
|
|
||||||
Generating classification report...
|
|
||||||
Saving model...
|
|
||||||
[MODEL MANAGER] Saving model to: model_train/model_mobilenet-lraspp_auto.joblib
|
|
||||||
[MODEL MANAGER] Saving metadata to: model_train/model_mobilenet-lraspp_auto_info.json
|
|
||||||
[MODEL MANAGER] Model saved successfully!
|
|
||||||
Training complete!
|
|
||||||
✅ Hoàn thành! Accuracy: 0.5669
|
|
||||||
@@ -1,15 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU HUẤN LUYỆN MÔ HÌNH RANDOM FOREST (GPU & CACHE)
|
|
||||||
Initializing FeatureExtractor (mode=extended)...
|
|
||||||
📦 Đang load cache: training_data_507bd2ba4ec0d3fe107839cbf73a7a7d.joblib...
|
|
||||||
✅ Loaded 632 samples từ cache!
|
|
||||||
⚡ Đã bỏ qua download위성 data (tiết kiệm thời gian)
|
|
||||||
[CACHE HIT] Using cached dataset with 632 samples
|
|
||||||
Training RANDOM_FOREST model...
|
|
||||||
Evaluating model...
|
|
||||||
Generating classification report...
|
|
||||||
Saving model...
|
|
||||||
[MODEL MANAGER] Saving model to: model_train/model_random_forest_auto.joblib
|
|
||||||
[MODEL MANAGER] Saving metadata to: model_train/model_random_forest_auto_info.json
|
|
||||||
[MODEL MANAGER] Model saved successfully!
|
|
||||||
Training complete!
|
|
||||||
✅ Hoàn thành! Accuracy: 0.6142
|
|
||||||
@@ -1,15 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU HUẤN LUYỆN MÔ HÌNH SVM (GPU & CACHE)
|
|
||||||
Initializing FeatureExtractor (mode=extended)...
|
|
||||||
📦 Đang load cache: training_data_507bd2ba4ec0d3fe107839cbf73a7a7d.joblib...
|
|
||||||
✅ Loaded 632 samples từ cache!
|
|
||||||
⚡ Đã bỏ qua download위성 data (tiết kiệm thời gian)
|
|
||||||
[CACHE HIT] Using cached dataset with 632 samples
|
|
||||||
Training SVM model...
|
|
||||||
Evaluating model...
|
|
||||||
Generating classification report...
|
|
||||||
Saving model...
|
|
||||||
[MODEL MANAGER] Saving model to: model_train/model_svm_auto.joblib
|
|
||||||
[MODEL MANAGER] Saving metadata to: model_train/model_svm_auto_info.json
|
|
||||||
[MODEL MANAGER] Model saved successfully!
|
|
||||||
Training complete!
|
|
||||||
✅ Hoàn thành! Accuracy: 0.5827
|
|
||||||
@@ -1,30 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU HUẤN LUYỆN MÔ HÌNH SWIN-UNET (GPU & CACHE)
|
|
||||||
Initializing FeatureExtractor (mode=extended)...
|
|
||||||
📦 Đang load cache: training_data_507bd2ba4ec0d3fe107839cbf73a7a7d.joblib...
|
|
||||||
✅ Loaded 632 samples từ cache!
|
|
||||||
⚡ Đã bỏ qua download위성 data (tiết kiệm thời gian)
|
|
||||||
[CACHE HIT] Using cached dataset with 632 samples
|
|
||||||
Training SWIN-UNET model...
|
|
||||||
Building Swin-UNet model on cuda...
|
|
||||||
[SWIN-UNET] Class distribution: [ 48 89 3 86 74 38 117 50]
|
|
||||||
[SWIN-UNET] Class weights: [0.37418982 0.20181024 5.98703525 0.20885013 0.24271772 0.47266082
|
|
||||||
0.15351377 0.35922223]
|
|
||||||
Training Swin-UNet model with PyTorch (with class weights)...
|
|
||||||
Swin-UNet Epoch 5/40, Train Loss: 1.4096, Val Loss: 1.2804, Val Acc: 48.03%, LR: 0.000293
|
|
||||||
[SWIN-UNET] Epoch 5/40 - Train Loss: 1.4096, Val Loss: 1.2804, Val Acc: 48.03%
|
|
||||||
Swin-UNet Epoch 10/40, Train Loss: 1.1129, Val Loss: 1.2746, Val Acc: 57.48%, LR: 0.000271
|
|
||||||
[SWIN-UNET] Epoch 10/40 - Train Loss: 1.1129, Val Loss: 1.2746, Val Acc: 57.48%
|
|
||||||
Swin-UNet Epoch 15/40, Train Loss: 1.0479, Val Loss: 1.3126, Val Acc: 50.39%, LR: 0.000238
|
|
||||||
[SWIN-UNET] Epoch 15/40 - Train Loss: 1.0479, Val Loss: 1.3126, Val Acc: 50.39%
|
|
||||||
Swin-UNet Epoch 20/40, Train Loss: 0.9548, Val Loss: 1.1118, Val Acc: 52.76%, LR: 0.000196
|
|
||||||
[SWIN-UNET] Epoch 20/40 - Train Loss: 0.9548, Val Loss: 1.1118, Val Acc: 52.76%
|
|
||||||
[SWIN-UNET] Early stopping at epoch 22 (best val loss: 1.1096)
|
|
||||||
Swin-UNet early stopped at epoch 22
|
|
||||||
Evaluating model...
|
|
||||||
Generating classification report...
|
|
||||||
Saving model...
|
|
||||||
[MODEL MANAGER] Saving model to: model_train/model_swin-unet_auto.joblib
|
|
||||||
[MODEL MANAGER] Saving metadata to: model_train/model_swin-unet_auto_info.json
|
|
||||||
[MODEL MANAGER] Model saved successfully!
|
|
||||||
Training complete!
|
|
||||||
✅ Hoàn thành! Accuracy: 0.5354
|
|
||||||
-94521
File diff suppressed because it is too large
Load Diff
@@ -1,38 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU PIPELINE 2D PATCH-BASED & CLOUD REMOVAL
|
|
||||||
Loading 2D patches from dataset_cache/training_data_2d.joblib...
|
|
||||||
Training 2D CNN with Data Augmentation...
|
|
||||||
Epoch 1/150 - Loss: 2.2272 - Test Acc: 0.0752 🌟
|
|
||||||
Epoch 2/150 - Loss: 2.0843 - Test Acc: 0.0796 🌟
|
|
||||||
Epoch 4/150 - Loss: 2.0350 - Test Acc: 0.1372 🌟
|
|
||||||
Epoch 8/150 - Loss: 2.0447 - Test Acc: 0.1637 🌟
|
|
||||||
Epoch 10/150 - Loss: 2.0397 - Test Acc: 0.1372
|
|
||||||
Epoch 12/150 - Loss: 2.0353 - Test Acc: 0.2168 🌟
|
|
||||||
Epoch 20/150 - Loss: 2.0139 - Test Acc: 0.1372
|
|
||||||
Traceback (most recent call last):
|
|
||||||
File "/home/x79/remote-sensing/train_land_2d_patch.py", line 298, in <module>
|
|
||||||
main()
|
|
||||||
File "/home/x79/remote-sensing/train_land_2d_patch.py", line 294, in main
|
|
||||||
train_2d_model(X, y)
|
|
||||||
File "/home/x79/remote-sensing/train_land_2d_patch.py", line 219, in train_2d_model
|
|
||||||
for batch_X, batch_y in train_loader:
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torch/utils/data/dataloader.py", line 725, in __next__
|
|
||||||
data = self._next_data()
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torch/utils/data/dataloader.py", line 785, in _next_data
|
|
||||||
data = self._dataset_fetcher.fetch(index) # may raise StopIteration
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torch/utils/data/_utils/fetch.py", line 54, in fetch
|
|
||||||
data = [self.dataset[idx] for idx in possibly_batched_index]
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torch/utils/data/_utils/fetch.py", line 54, in <listcomp>
|
|
||||||
data = [self.dataset[idx] for idx in possibly_batched_index]
|
|
||||||
File "/home/x79/remote-sensing/train_land_2d_patch.py", line 192, in __getitem__
|
|
||||||
x = transform(x)
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torchvision/transforms/transforms.py", line 95, in __call__
|
|
||||||
img = t(img)
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torch/nn/modules/module.py", line 1778, in _wrapped_call_impl
|
|
||||||
return self._call_impl(*args, **kwargs)
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torch/nn/modules/module.py", line 1789, in _call_impl
|
|
||||||
return forward_call(*args, **kwargs)
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torchvision/transforms/transforms.py", line 752, in forward
|
|
||||||
return F.vflip(img)
|
|
||||||
File "/home/x79/miniconda3/envs/env_01/lib/python3.10/site-packages/torchvision/transforms/functional.py", line 757, in vflip
|
|
||||||
def vflip(img: Tensor) -> Tensor:
|
|
||||||
KeyboardInterrupt
|
|
||||||
File diff suppressed because it is too large
Load Diff
File diff suppressed because one or more lines are too long
@@ -1,39 +0,0 @@
|
|||||||
🚀 V4: TÍCH HỢP RADAR SENTINEL-1 (32-CHANNELS FUSION)
|
|
||||||
============================================================
|
|
||||||
Clean FUSION data: (252, 32, 16, 16), 7 classes, [31, 38, 32, 49, 23, 75, 4]
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
32-CHANNELS FUSION CNN
|
|
||||||
============================================================
|
|
||||||
Ep 1 Fusion-Acc=0.1961 🌟
|
|
||||||
Ep 3 Fusion-Acc=0.2745 🌟
|
|
||||||
Ep 4 Fusion-Acc=0.4706 🌟
|
|
||||||
Ep 5 Fusion-Acc=0.5490 🌟
|
|
||||||
Ep 6 Fusion-Acc=0.7059 🌟
|
|
||||||
Ep 7 Fusion-Acc=0.7451 🌟
|
|
||||||
Ep 8 Fusion-Acc=0.8431 🌟
|
|
||||||
Ep 15 Fusion-Acc=0.8627 🌟
|
|
||||||
|
|
||||||
✅ CNN Fusion best: 0.8627
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
HYBRID FUSION: CNN embed + S1/S2 Rich features + XGBoost
|
|
||||||
============================================================
|
|
||||||
Extracted 2182 fusion features per sample
|
|
||||||
Final Feature Vector: (252, 2694)
|
|
||||||
✅ Hybrid Fusion Acc: 0.8627
|
|
||||||
Fold 1: 0.9412
|
|
||||||
Fold 2: 0.9020
|
|
||||||
Fold 3: 0.9200
|
|
||||||
Fold 4: 0.8800
|
|
||||||
Fold 5: 0.9000
|
|
||||||
✅ CV Mean: 0.9086 ± 0.0206
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
📊 FINAL RESULTS V4 (WITH RADAR)
|
|
||||||
============================================================
|
|
||||||
✅ Hybrid Fusion CV: 0.9086
|
|
||||||
📈 CNN Fusion (32ch): 0.8627
|
|
||||||
📈 Hybrid Fusion (CNN+XGB): 0.8627
|
|
||||||
|
|
||||||
🏆 BEST: 0.9086
|
|
||||||
@@ -1,18 +0,0 @@
|
|||||||
|
|
||||||
============================================================
|
|
||||||
HYBRID FUSION ENSEMBLE: CNN embed + S1/S2 Rich features + XGB/LGBM/ETC
|
|
||||||
============================================================
|
|
||||||
Final Feature Vector: (443, 2694)
|
|
||||||
Fold 1: 0.8876
|
|
||||||
Fold 2: 0.9438
|
|
||||||
Fold 3: 0.9438
|
|
||||||
Fold 4: 0.9659
|
|
||||||
Fold 5: 0.9432
|
|
||||||
✅ Ensemble CV Mean: 0.9369 ± 0.0261
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
📊 FINAL RESULTS V5 (ENSEMBLE + RADAR)
|
|
||||||
============================================================
|
|
||||||
✅ Hybrid Fusion Ensemble CV: 0.9369
|
|
||||||
|
|
||||||
🏆 BEST: 0.9369
|
|
||||||
@@ -1,21 +0,0 @@
|
|||||||
🚀 V6: EXHAUSTIVE HYPERPARAMETER TUNING
|
|
||||||
============================================================
|
|
||||||
Data: (443, 32, 16, 16), 7 classes, dist=[65, 55, 48, 72, 74, 124, 5]
|
|
||||||
|
|
||||||
--- Training Multi-Seed CNN Ensemble ---
|
|
||||||
Seed 42: CNN Acc = 0.8989
|
|
||||||
Seed 123: CNN Acc = 0.8652
|
|
||||||
Seed 777: CNN Acc = 0.8876
|
|
||||||
Multi-seed CNN embedding: (443, 1536)
|
|
||||||
Total features: (443, 3940)
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
🔬 EXHAUSTIVE HYPERPARAMETER SEARCH
|
|
||||||
============================================================
|
|
||||||
🏆 XGB-deep: 0.9526 ± 0.0110 (folds: ['0.955', '0.933', '0.966', '0.955', '0.955'])
|
|
||||||
✅ XGB-shallow: 0.9436 ± 0.0173 (folds: ['0.944', '0.910', '0.955', '0.955', '0.955'])
|
|
||||||
🏆 XGB-balanced: 0.9504 ± 0.0113 (folds: ['0.944', '0.933', '0.955', '0.955', '0.966'])
|
|
||||||
✅ LGBM-tuned: 0.9458 ± 0.0149 (folds: ['0.955', '0.921', '0.966', '0.943', '0.943'])
|
|
||||||
✅ LGBM-conservative: 0.9481 ± 0.0152 (folds: ['0.944', '0.921', '0.966', '0.955', '0.955'])
|
|
||||||
🏆 ETC-deep: 0.9572 ± 0.0082 (folds: ['0.944', '0.955', '0.955', '0.966', '0.966'])
|
|
||||||
🏆 RF-tuned: 0.9549 ± 0.0099 (folds: ['0.944', '0.955', '0.944', '0.966', '0.966'])
|
|
||||||
@@ -1,157 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU TÌM KIẾM SIÊU THAM SỐ CHO SWIN-UNET
|
|
||||||
Loading data from dataset_cache/training_data_507bd2ba4ec0d3fe107839cbf73a7a7d.joblib...
|
|
||||||
Using device: cuda
|
|
||||||
|
|
||||||
[1/48] Training with params: {'embed_dim': 64, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.5984
|
|
||||||
🌟 NEW BEST ACCURACY: 0.5984
|
|
||||||
|
|
||||||
[2/48] Training with params: {'embed_dim': 64, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6063
|
|
||||||
🌟 NEW BEST ACCURACY: 0.6063
|
|
||||||
|
|
||||||
[3/48] Training with params: {'embed_dim': 64, 'lr': 0.001, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6142
|
|
||||||
🌟 NEW BEST ACCURACY: 0.6142
|
|
||||||
|
|
||||||
[4/48] Training with params: {'embed_dim': 64, 'lr': 0.001, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6772
|
|
||||||
🌟 NEW BEST ACCURACY: 0.6772
|
|
||||||
|
|
||||||
[5/48] Training with params: {'embed_dim': 64, 'lr': 0.0005, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.5827
|
|
||||||
|
|
||||||
[6/48] Training with params: {'embed_dim': 64, 'lr': 0.0005, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.5984
|
|
||||||
|
|
||||||
[7/48] Training with params: {'embed_dim': 64, 'lr': 0.0005, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.5669
|
|
||||||
|
|
||||||
[8/48] Training with params: {'embed_dim': 64, 'lr': 0.0005, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6142
|
|
||||||
|
|
||||||
[9/48] Training with params: {'embed_dim': 64, 'lr': 0.0001, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6142
|
|
||||||
|
|
||||||
[10/48] Training with params: {'embed_dim': 64, 'lr': 0.0001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.5118
|
|
||||||
|
|
||||||
[11/48] Training with params: {'embed_dim': 64, 'lr': 0.0001, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.5591
|
|
||||||
|
|
||||||
[12/48] Training with params: {'embed_dim': 64, 'lr': 0.0001, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6063
|
|
||||||
|
|
||||||
[13/48] Training with params: {'embed_dim': 128, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.7008
|
|
||||||
🌟 NEW BEST ACCURACY: 0.7008
|
|
||||||
|
|
||||||
[14/48] Training with params: {'embed_dim': 128, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6142
|
|
||||||
|
|
||||||
[15/48] Training with params: {'embed_dim': 128, 'lr': 0.001, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6772
|
|
||||||
|
|
||||||
[16/48] Training with params: {'embed_dim': 128, 'lr': 0.001, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6063
|
|
||||||
|
|
||||||
[17/48] Training with params: {'embed_dim': 128, 'lr': 0.0005, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.5906
|
|
||||||
|
|
||||||
[18/48] Training with params: {'embed_dim': 128, 'lr': 0.0005, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6457
|
|
||||||
|
|
||||||
[19/48] Training with params: {'embed_dim': 128, 'lr': 0.0005, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.5906
|
|
||||||
|
|
||||||
[20/48] Training with params: {'embed_dim': 128, 'lr': 0.0005, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.5906
|
|
||||||
|
|
||||||
[21/48] Training with params: {'embed_dim': 128, 'lr': 0.0001, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6220
|
|
||||||
|
|
||||||
[22/48] Training with params: {'embed_dim': 128, 'lr': 0.0001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6457
|
|
||||||
|
|
||||||
[23/48] Training with params: {'embed_dim': 128, 'lr': 0.0001, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.5984
|
|
||||||
|
|
||||||
[24/48] Training with params: {'embed_dim': 128, 'lr': 0.0001, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6063
|
|
||||||
|
|
||||||
[25/48] Training with params: {'embed_dim': 256, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.7087
|
|
||||||
🌟 NEW BEST ACCURACY: 0.7087
|
|
||||||
|
|
||||||
[26/48] Training with params: {'embed_dim': 256, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.7323
|
|
||||||
🌟 NEW BEST ACCURACY: 0.7323
|
|
||||||
|
|
||||||
[27/48] Training with params: {'embed_dim': 256, 'lr': 0.001, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6457
|
|
||||||
|
|
||||||
[28/48] Training with params: {'embed_dim': 256, 'lr': 0.001, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6142
|
|
||||||
|
|
||||||
[29/48] Training with params: {'embed_dim': 256, 'lr': 0.0005, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.5984
|
|
||||||
|
|
||||||
[30/48] Training with params: {'embed_dim': 256, 'lr': 0.0005, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.5906
|
|
||||||
|
|
||||||
[31/48] Training with params: {'embed_dim': 256, 'lr': 0.0005, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6142
|
|
||||||
|
|
||||||
[32/48] Training with params: {'embed_dim': 256, 'lr': 0.0005, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.7008
|
|
||||||
|
|
||||||
[33/48] Training with params: {'embed_dim': 256, 'lr': 0.0001, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6299
|
|
||||||
|
|
||||||
[34/48] Training with params: {'embed_dim': 256, 'lr': 0.0001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6299
|
|
||||||
|
|
||||||
[35/48] Training with params: {'embed_dim': 256, 'lr': 0.0001, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6378
|
|
||||||
|
|
||||||
[36/48] Training with params: {'embed_dim': 256, 'lr': 0.0001, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6457
|
|
||||||
|
|
||||||
[37/48] Training with params: {'embed_dim': 512, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6142
|
|
||||||
|
|
||||||
[38/48] Training with params: {'embed_dim': 512, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6457
|
|
||||||
|
|
||||||
[39/48] Training with params: {'embed_dim': 512, 'lr': 0.001, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6378
|
|
||||||
|
|
||||||
[40/48] Training with params: {'embed_dim': 512, 'lr': 0.001, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6299
|
|
||||||
|
|
||||||
[41/48] Training with params: {'embed_dim': 512, 'lr': 0.0005, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6378
|
|
||||||
|
|
||||||
[42/48] Training with params: {'embed_dim': 512, 'lr': 0.0005, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6220
|
|
||||||
|
|
||||||
[43/48] Training with params: {'embed_dim': 512, 'lr': 0.0005, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6142
|
|
||||||
|
|
||||||
[44/48] Training with params: {'embed_dim': 512, 'lr': 0.0005, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.7008
|
|
||||||
|
|
||||||
[45/48] Training with params: {'embed_dim': 512, 'lr': 0.0001, 'weight_decay': 0.01, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6457
|
|
||||||
|
|
||||||
[46/48] Training with params: {'embed_dim': 512, 'lr': 0.0001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.7323
|
|
||||||
|
|
||||||
[47/48] Training with params: {'embed_dim': 512, 'lr': 0.0001, 'weight_decay': 0.001, 'epochs': 200}
|
|
||||||
Test Accuracy: 0.6693
|
|
||||||
|
|
||||||
[48/48] Training with params: {'embed_dim': 512, 'lr': 0.0001, 'weight_decay': 0.001, 'epochs': 500}
|
|
||||||
Test Accuracy: 0.6457
|
|
||||||
|
|
||||||
✅ Đã lưu mô hình tốt nhất (Acc: 0.7323) vào land_classification_model/model_swin-unet_optimized_95.joblib
|
|
||||||
Cấu hình tốt nhất: {'embed_dim': 256, 'lr': 0.001, 'weight_decay': 0.01, 'epochs': 500}
|
|
||||||
@@ -1,117 +0,0 @@
|
|||||||
🚀 CHIẾN LƯỢC TOÀN DIỆN ĐẠT >95% ACCURACY
|
|
||||||
============================================================
|
|
||||||
Loaded data: X=(706, 24, 16, 16), y=(706,)
|
|
||||||
Labels unique: [-1 0 1 2 3 4 5 6]
|
|
||||||
After cleanup: X=(652, 24, 16, 16), y=(652,) (removed 54 bad samples)
|
|
||||||
Remapped labels: [0 1 2 3 4 5 6]
|
|
||||||
Class 0: 65 samples
|
|
||||||
Class 1: 52 samples
|
|
||||||
Class 2: 48 samples
|
|
||||||
Class 3: 72 samples
|
|
||||||
Class 4: 108 samples
|
|
||||||
Class 5: 219 samples
|
|
||||||
Class 6: 88 samples
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
STRATEGY 5: Flat pixel features + XGBoost (sanity check)
|
|
||||||
============================================================
|
|
||||||
Flat features: (652, 6144)
|
|
||||||
✅ Flat XGBoost acc: 0.7328
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
STRATEGY 1: Lightweight CNN (no upsampling)
|
|
||||||
============================================================
|
|
||||||
Device: cuda
|
|
||||||
Epoch 1/300 Loss=1.8379 Acc=0.0763 🌟
|
|
||||||
Epoch 2/300 Loss=1.5825 Acc=0.2824 🌟
|
|
||||||
Epoch 3/300 Loss=1.4412 Acc=0.5649 🌟
|
|
||||||
Epoch 4/300 Loss=1.3274 Acc=0.6336 🌟
|
|
||||||
Epoch 5/300 Loss=1.2893 Acc=0.6870 🌟
|
|
||||||
Epoch 8/300 Loss=1.2140 Acc=0.7099 🌟
|
|
||||||
Epoch 11/300 Loss=1.2224 Acc=0.7252 🌟
|
|
||||||
Epoch 14/300 Loss=1.1087 Acc=0.7328 🌟
|
|
||||||
Epoch 16/300 Loss=1.1081 Acc=0.7710 🌟
|
|
||||||
Epoch 18/300 Loss=1.1847 Acc=0.7939 🌟
|
|
||||||
Epoch 20/300 Loss=1.0316 Acc=0.7786 (patience=2)
|
|
||||||
Epoch 24/300 Loss=1.0465 Acc=0.8092 🌟
|
|
||||||
Epoch 26/300 Loss=0.9583 Acc=0.8397 🌟
|
|
||||||
Epoch 40/300 Loss=0.9361 Acc=0.8626 🌟
|
|
||||||
Epoch 60/300 Loss=0.9807 Acc=0.8092 (patience=20)
|
|
||||||
Epoch 80/300 Loss=0.9459 Acc=0.8015 (patience=40)
|
|
||||||
Epoch 100/300 Loss=0.8443 Acc=0.7939 (patience=60)
|
|
||||||
Early stop at epoch 100
|
|
||||||
✅ LightCNN best acc: 0.8626
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
STRATEGY 2: Hybrid CNN embeddings + XGBoost
|
|
||||||
============================================================
|
|
||||||
CNN embeddings: (652, 256)
|
|
||||||
Extracted 316 rich features per sample
|
|
||||||
Combined features: (652, 572)
|
|
||||||
✅ Hybrid XGBoost acc: 0.8244
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
STRATEGY 3: Rich Features + Stacking Ensemble
|
|
||||||
============================================================
|
|
||||||
Extracted 316 rich features per sample
|
|
||||||
XGBoost: 0.7939
|
|
||||||
LightGBM: 0.7786
|
|
||||||
ExtraTrees: 0.7863
|
|
||||||
RandomForest: 0.7710
|
|
||||||
GBM: 0.7710
|
|
||||||
[06:55:44] WARNING: /__w/xgboost/xgboost/src/learner.cc:782:
|
|
||||||
Parameters: { "use_label_encoder" } are not used.
|
|
||||||
[06:55:44] WARNING: /__w/xgboost/xgboost/src/learner.cc:782:
|
|
||||||
Parameters: { "use_label_encoder" } are not used.
|
|
||||||
|
|
||||||
|
|
||||||
[06:55:44] WARNING: /__w/xgboost/xgboost/src/learner.cc:782:
|
|
||||||
Parameters: { "use_label_encoder" } are not used.
|
|
||||||
[06:55:44] WARNING: /__w/xgboost/xgboost/src/learner.cc:782:
|
|
||||||
Parameters: { "use_label_encoder" } are not used.
|
|
||||||
|
|
||||||
|
|
||||||
[06:55:44] WARNING: /__w/xgboost/xgboost/src/learner.cc:782:
|
|
||||||
Parameters: { "use_label_encoder" } are not used.
|
|
||||||
|
|
||||||
[06:56:17] WARNING: /__w/xgboost/xgboost/src/common/error_msg.cc:62: Falling back to prediction using DMatrix due to mismatched devices. This might lead to higher memory usage and slower performance. XGBoost is running on: cuda:0, while the input data is on: cpu.
|
|
||||||
Potential solutions:
|
|
||||||
- Use a data structure that matches the device ordinal in the booster.
|
|
||||||
- Set the device for booster before call to inplace_predict.
|
|
||||||
|
|
||||||
This warning will only be shown once.
|
|
||||||
|
|
||||||
Stacking Ensemble: 0.7710
|
|
||||||
Voting Ensemble: 0.7786
|
|
||||||
✅ Best ensemble: XGBoost = 0.7939
|
|
||||||
Extracted 316 rich features per sample
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
STRATEGY 4: 5-Fold Stratified Cross-Validation
|
|
||||||
============================================================
|
|
||||||
Fold 1: 0.7939
|
|
||||||
Fold 2: 0.8244
|
|
||||||
Fold 3: 0.7769
|
|
||||||
Fold 4: 0.8154
|
|
||||||
Fold 5: 0.7615
|
|
||||||
✅ CV Mean: 0.7944 ± 0.0234
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
📊 TỔNG KẾT KẾT QUẢ
|
|
||||||
============================================================
|
|
||||||
📈 LightCNN: 0.8626
|
|
||||||
📈 Hybrid CNN+XGBoost: 0.8244
|
|
||||||
📈 CV Mean (XGBoost rich): 0.7944
|
|
||||||
📈 Ensemble XGBoost: 0.7939
|
|
||||||
📈 Ensemble ExtraTrees: 0.7863
|
|
||||||
📈 Ensemble LightGBM: 0.7786
|
|
||||||
📈 Ensemble Voting: 0.7786
|
|
||||||
📈 Ensemble RandomForest: 0.7710
|
|
||||||
📈 Ensemble GBM: 0.7710
|
|
||||||
📈 Ensemble Stacking: 0.7710
|
|
||||||
📈 Flat XGBoost (baseline): 0.7328
|
|
||||||
|
|
||||||
🏆 BEST: LightCNN = 0.8626
|
|
||||||
|
|
||||||
✅ Kết quả đã được lưu vào model_train/ultimate_results.json
|
|
||||||
⚠️ Chưa đạt 95%. Best = 0.8626. Cần thêm dữ liệu hoặc feature engineering.
|
|
||||||
@@ -1,69 +0,0 @@
|
|||||||
🚀 CHIẾN LƯỢC V2: TOÀN DIỆN ĐẠT >95%
|
|
||||||
============================================================
|
|
||||||
Clean data: (652, 24, 16, 16), 7 classes
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
CNN + TTA (Test-Time Augmentation)
|
|
||||||
============================================================
|
|
||||||
Ep 1 Loss=1.6500 TTA-Acc=0.1450 🌟
|
|
||||||
Ep 2 Loss=1.4633 TTA-Acc=0.3511 🌟
|
|
||||||
Ep 3 Loss=1.3316 TTA-Acc=0.6336 🌟
|
|
||||||
Ep 4 Loss=1.2101 TTA-Acc=0.7252 🌟
|
|
||||||
Ep 6 Loss=1.2056 TTA-Acc=0.8015 🌟
|
|
||||||
Ep 13 Loss=1.1216 TTA-Acc=0.8244 🌟
|
|
||||||
Ep 23 Loss=0.9355 TTA-Acc=0.8321 🌟
|
|
||||||
Ep 30 Loss=0.9346 TTA-Acc=0.8092 (pat=7)
|
|
||||||
Ep 60 Loss=0.9658 TTA-Acc=0.7634 (pat=37)
|
|
||||||
Ep 69 Loss=0.8988 TTA-Acc=0.8397 🌟
|
|
||||||
Ep 74 Loss=0.9159 TTA-Acc=0.8550 🌟
|
|
||||||
Ep 90 Loss=0.8761 TTA-Acc=0.8168 (pat=16)
|
|
||||||
Ep 120 Loss=0.9473 TTA-Acc=0.8092 (pat=46)
|
|
||||||
Ep 139 Loss=0.9317 TTA-Acc=0.8626 🌟
|
|
||||||
Ep 150 Loss=0.7716 TTA-Acc=0.8397 (pat=11)
|
|
||||||
Ep 175 Loss=0.7432 TTA-Acc=0.8702 🌟
|
|
||||||
Ep 180 Loss=0.8290 TTA-Acc=0.8702 (pat=5)
|
|
||||||
Ep 210 Loss=0.8060 TTA-Acc=0.8626 (pat=35)
|
|
||||||
Ep 240 Loss=0.6980 TTA-Acc=0.8473 (pat=65)
|
|
||||||
Early stop ep 255
|
|
||||||
✅ CNN+TTA best: 0.8702
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
RICH FEATURES V2 + ENSEMBLE
|
|
||||||
============================================================
|
|
||||||
Extracted 1713 features per sample
|
|
||||||
XGB: 0.7557
|
|
||||||
LGBM: 0.7634
|
|
||||||
ET: 0.7481
|
|
||||||
RF: 0.7557
|
|
||||||
Voting: 0.7634
|
|
||||||
|
|
||||||
5-Fold CV:
|
|
||||||
Fold 1: 0.7557
|
|
||||||
Fold 2: 0.7939
|
|
||||||
Fold 3: 0.7769
|
|
||||||
Fold 4: 0.8308
|
|
||||||
Fold 5: 0.8000
|
|
||||||
CV: 0.7915 ± 0.0249
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
HYBRID V2: CNN embed + Rich features + XGBoost
|
|
||||||
============================================================
|
|
||||||
Extracted 1713 features per sample
|
|
||||||
Combined: (652, 2097)
|
|
||||||
✅ Hybrid V2: 0.8244
|
|
||||||
CV: 0.9142 ± 0.0194
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
📊 KẾT QUẢ TỔNG HỢP V2
|
|
||||||
============================================================
|
|
||||||
✅ Hybrid CV: 0.9142
|
|
||||||
📈 CNN+TTA: 0.8702
|
|
||||||
📈 Hybrid V2: 0.8244
|
|
||||||
📈 Ens CV: 0.7915
|
|
||||||
📈 Ens_LGBM: 0.7634
|
|
||||||
📈 Ens_Vote: 0.7634
|
|
||||||
📈 Ens_XGB: 0.7557
|
|
||||||
📈 Ens_RF: 0.7557
|
|
||||||
📈 Ens_ET: 0.7481
|
|
||||||
|
|
||||||
🏆 BEST: Hybrid CV = 0.9142
|
|
||||||
@@ -1,27 +0,0 @@
|
|||||||
🚀 V3: MULTI-SEED ENSEMBLE + T0-ONLY + SELF-TRAINING
|
|
||||||
============================================================
|
|
||||||
Clean: (652, 24, 16, 16), 7 classes, [65, 52, 48, 72, 108, 219, 88]
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
MULTI-SEED CNN ENSEMBLE (10 models)
|
|
||||||
============================================================
|
|
||||||
Seed 0: 0.8473
|
|
||||||
Seed 1: 0.8702
|
|
||||||
Seed 2: 0.8550
|
|
||||||
Seed 3: 0.8550
|
|
||||||
Seed 4: 0.8702
|
|
||||||
Seed 5: 0.8626
|
|
||||||
Seed 6: 0.8702
|
|
||||||
Seed 7: 0.8702
|
|
||||||
Seed 8: 0.8702
|
|
||||||
Seed 9: 0.8702
|
|
||||||
✅ 10-Model Ensemble TTA: 0.8473
|
|
||||||
|
|
||||||
============================================================
|
|
||||||
TIMESTEP-0-ONLY XGBoost (cleanest data)
|
|
||||||
============================================================
|
|
||||||
T0 valid: 539/652
|
|
||||||
Features: (539, 1638)
|
|
||||||
XGB t0: 0.6852
|
|
||||||
LGBM t0: 0.6296
|
|
||||||
ET t0: 0.6852
|
|
||||||
@@ -1,15 +0,0 @@
|
|||||||
🚀 BẮT ĐẦU HUẤN LUYỆN MÔ HÌNH XGBOOST (GPU & CACHE)
|
|
||||||
Initializing FeatureExtractor (mode=extended)...
|
|
||||||
📦 Đang load cache: training_data_507bd2ba4ec0d3fe107839cbf73a7a7d.joblib...
|
|
||||||
✅ Loaded 632 samples từ cache!
|
|
||||||
⚡ Đã bỏ qua download위성 data (tiết kiệm thời gian)
|
|
||||||
[CACHE HIT] Using cached dataset with 632 samples
|
|
||||||
Training XGBOOST model...
|
|
||||||
Evaluating model...
|
|
||||||
Generating classification report...
|
|
||||||
Saving model...
|
|
||||||
[MODEL MANAGER] Saving model to: model_train/model_xgboost_auto.joblib
|
|
||||||
[MODEL MANAGER] Saving metadata to: model_train/model_xgboost_auto_info.json
|
|
||||||
[MODEL MANAGER] Model saved successfully!
|
|
||||||
Training complete!
|
|
||||||
✅ Hoàn thành! Accuracy: 0.6378
|
|
||||||
File diff suppressed because one or more lines are too long
@@ -1,231 +0,0 @@
|
|||||||
import re
|
|
||||||
import os
|
|
||||||
|
|
||||||
with open('/home/x79/remote-sensing/new_import_ODC.py', 'r', encoding='utf-8') as f:
|
|
||||||
content = f.read()
|
|
||||||
|
|
||||||
# 1. load_data
|
|
||||||
load_data_replacement = """def load_data(dc, date_range, longtitude_range, latitude_range):
|
|
||||||
import os, hashlib
|
|
||||||
cache_dir = "dataset_cache"
|
|
||||||
os.makedirs(cache_dir, exist_ok=True)
|
|
||||||
key_str = f"s2_{date_range}_{longtitude_range}_{latitude_range}"
|
|
||||||
cache_key = hashlib.md5(key_str.encode()).hexdigest() + ".nc"
|
|
||||||
cache_path = os.path.join(cache_dir, cache_key)
|
|
||||||
|
|
||||||
if os.path.exists(cache_path):
|
|
||||||
print(f"✅ Loading cached S2 data from {cache_path}")
|
|
||||||
return xr.open_dataset(cache_path)
|
|
||||||
|
|
||||||
product = 's2_l2a'
|
|
||||||
bbox = [longtitude_range[0], latitude_range[0], longtitude_range[1], latitude_range[1]]
|
|
||||||
|
|
||||||
import pystac_client
|
|
||||||
import planetary_computer
|
|
||||||
import odc.stac
|
|
||||||
|
|
||||||
catalog = pystac_client.Client.open(
|
|
||||||
"https://planetarycomputer.microsoft.com/api/stac/v1",
|
|
||||||
modifier=planetary_computer.sign_inplace,
|
|
||||||
)
|
|
||||||
search = catalog.search(
|
|
||||||
collections=["sentinel-2-l2a"],
|
|
||||||
bbox=bbox,
|
|
||||||
datetime=f"{date_range[0]}/{date_range[1]}",
|
|
||||||
)
|
|
||||||
items = list(search.items())
|
|
||||||
|
|
||||||
data = odc.stac.load(
|
|
||||||
items,
|
|
||||||
bands=["red", "nir", "SCL"],
|
|
||||||
bbox=bbox,
|
|
||||||
crs="EPSG:32648",
|
|
||||||
resolution=RESOLUTION,
|
|
||||||
chunks={"x": 2048, "y": 2048, "time": 1},
|
|
||||||
groupby="solar_day"
|
|
||||||
)
|
|
||||||
if "SCL" in data.data_vars:
|
|
||||||
data = data.rename({"SCL": "scl"})
|
|
||||||
|
|
||||||
print(f"💾 Caching S2 data to {cache_path}")
|
|
||||||
data = data.compute()
|
|
||||||
data.to_netcdf(cache_path)
|
|
||||||
|
|
||||||
return data"""
|
|
||||||
content = re.sub(r'def load_data\(dc, date_range, longtitude_range, latitude_range\):.*?return data', load_data_replacement, content, flags=re.DOTALL)
|
|
||||||
|
|
||||||
|
|
||||||
# 2. load_sen1
|
|
||||||
load_sen1_replacement = """def load_sen1(bbox, time_range):
|
|
||||||
import os, hashlib
|
|
||||||
cache_dir = "dataset_cache"
|
|
||||||
os.makedirs(cache_dir, exist_ok=True)
|
|
||||||
key_str = f"s1_vh_vv_{bbox}_{time_range}"
|
|
||||||
cache_key_vh = hashlib.md5((key_str + "vh").encode()).hexdigest() + ".nc"
|
|
||||||
cache_key_vv = hashlib.md5((key_str + "vv").encode()).hexdigest() + ".nc"
|
|
||||||
cache_path_vh = os.path.join(cache_dir, cache_key_vh)
|
|
||||||
cache_path_vv = os.path.join(cache_dir, cache_key_vv)
|
|
||||||
|
|
||||||
if os.path.exists(cache_path_vh) and os.path.exists(cache_path_vv):
|
|
||||||
print(f"✅ Loading cached S1 data from {cache_path_vh} and {cache_path_vv}")
|
|
||||||
return xr.open_dataarray(cache_path_vh), xr.open_dataarray(cache_path_vv)
|
|
||||||
|
|
||||||
import pystac_client
|
|
||||||
import planetary_computer
|
|
||||||
import odc.stac
|
|
||||||
|
|
||||||
# Kết nối STAC Client
|
|
||||||
catalog = pystac_client.Client.open(
|
|
||||||
"https://planetarycomputer.microsoft.com/api/stac/v1",
|
|
||||||
modifier=planetary_computer.sign_inplace,
|
|
||||||
)
|
|
||||||
|
|
||||||
# Tìm kiếm Items
|
|
||||||
search = catalog.search(
|
|
||||||
collections=["sentinel-1-rtc"],
|
|
||||||
bbox=bbox,
|
|
||||||
datetime=time_range,
|
|
||||||
)
|
|
||||||
items = list(search.items())
|
|
||||||
|
|
||||||
# Tải dữ liệu thành xarray Dataset
|
|
||||||
ds_s1 = odc.stac.load(
|
|
||||||
items,
|
|
||||||
bands=["vv", "vh"],
|
|
||||||
bbox=bbox,
|
|
||||||
crs="EPSG:32648",
|
|
||||||
resolution=RESOLUTION,
|
|
||||||
chunks={"x": 2048, "y": 2048, "time": 1}
|
|
||||||
)
|
|
||||||
|
|
||||||
# Tính giá trị trung vị theo thời gian
|
|
||||||
ds_median = ds_s1.median(dim="time").compute()
|
|
||||||
vv = ds_median["vv"]
|
|
||||||
vh = ds_median["vh"]
|
|
||||||
|
|
||||||
# Thêm chiều 'band' để giống hệt rioxarray
|
|
||||||
vv = vv.expand_dims(dim="band")
|
|
||||||
vh = vh.expand_dims(dim="band")
|
|
||||||
|
|
||||||
# Phục hồi metadata về toạ độ
|
|
||||||
vv = vv.rio.write_crs("EPSG:32648")
|
|
||||||
vh = vh.rio.write_crs("EPSG:32648")
|
|
||||||
|
|
||||||
print(f"💾 Caching S1 data to {cache_dir}")
|
|
||||||
vh.to_netcdf(cache_path_vh)
|
|
||||||
vv.to_netcdf(cache_path_vv)
|
|
||||||
|
|
||||||
return vh, vv"""
|
|
||||||
content = re.sub(r'def load_sen1\(bbox, time_range\):.*?return vh, vv', load_sen1_replacement, content, flags=re.DOTALL)
|
|
||||||
|
|
||||||
|
|
||||||
# 3. load_data_sen1
|
|
||||||
load_data_sen1_replacement = """def load_data_sen1(dc, date_range, coordinates):
|
|
||||||
import os, hashlib
|
|
||||||
longtitude_range, latitude_range = coordinates
|
|
||||||
bbox = [longtitude_range[0], latitude_range[0], longtitude_range[1], latitude_range[1]]
|
|
||||||
|
|
||||||
cache_dir = "dataset_cache"
|
|
||||||
os.makedirs(cache_dir, exist_ok=True)
|
|
||||||
key_str = f"data_sen1_{date_range}_{bbox}"
|
|
||||||
cache_key_vh = hashlib.md5((key_str + "vh").encode()).hexdigest() + ".nc"
|
|
||||||
cache_key_vv = hashlib.md5((key_str + "vv").encode()).hexdigest() + ".nc"
|
|
||||||
cache_path_vh = os.path.join(cache_dir, cache_key_vh)
|
|
||||||
cache_path_vv = os.path.join(cache_dir, cache_key_vv)
|
|
||||||
|
|
||||||
if os.path.exists(cache_path_vh) and os.path.exists(cache_path_vv):
|
|
||||||
print(f"✅ Loading cached S1 (coord) data")
|
|
||||||
return xr.open_dataarray(cache_path_vh), xr.open_dataarray(cache_path_vv)
|
|
||||||
|
|
||||||
import pystac_client
|
|
||||||
import planetary_computer
|
|
||||||
import odc.stac
|
|
||||||
|
|
||||||
catalog = pystac_client.Client.open(
|
|
||||||
"https://planetarycomputer.microsoft.com/api/stac/v1",
|
|
||||||
modifier=planetary_computer.sign_inplace,
|
|
||||||
)
|
|
||||||
search = catalog.search(
|
|
||||||
collections=["sentinel-1-rtc"],
|
|
||||||
bbox=bbox,
|
|
||||||
datetime=f"{date_range[0]}/{date_range[1]}",
|
|
||||||
)
|
|
||||||
items = list(search.items())
|
|
||||||
|
|
||||||
data_sen1 = odc.stac.load(
|
|
||||||
items,
|
|
||||||
bands=["vv", "vh"],
|
|
||||||
bbox=bbox,
|
|
||||||
crs="EPSG:32648",
|
|
||||||
resolution=RESOLUTION,
|
|
||||||
chunks={"x": 2048, "y": 2048, "time": 1},
|
|
||||||
groupby="solar_day"
|
|
||||||
)
|
|
||||||
|
|
||||||
data_sen1 = data_sen1.compute()
|
|
||||||
dsvh = data_sen1.vh
|
|
||||||
dsvv = data_sen1.vv
|
|
||||||
|
|
||||||
print(f"💾 Caching S1 (coord) data")
|
|
||||||
dsvh.to_netcdf(cache_path_vh)
|
|
||||||
dsvv.to_netcdf(cache_path_vv)
|
|
||||||
|
|
||||||
return dsvh, dsvv"""
|
|
||||||
content = re.sub(r'def load_data_sen1\(dc, date_range, coordinates\):.*?return dsvh, dsvv', load_data_sen1_replacement, content, flags=re.DOTALL)
|
|
||||||
|
|
||||||
|
|
||||||
# 4. load_data_sen2
|
|
||||||
load_data_sen2_replacement = """def load_data_sen2(dc, date_range, coordinates):
|
|
||||||
import os, hashlib
|
|
||||||
longtitude_range, latitude_range = coordinates
|
|
||||||
bbox = [longtitude_range[0], latitude_range[0], longtitude_range[1], latitude_range[1]]
|
|
||||||
|
|
||||||
cache_dir = "dataset_cache"
|
|
||||||
os.makedirs(cache_dir, exist_ok=True)
|
|
||||||
key_str = f"data_sen2_{date_range}_{bbox}"
|
|
||||||
cache_key = hashlib.md5(key_str.encode()).hexdigest() + ".nc"
|
|
||||||
cache_path = os.path.join(cache_dir, cache_key)
|
|
||||||
|
|
||||||
if os.path.exists(cache_path):
|
|
||||||
print(f"✅ Loading cached S2 (coord) data from {cache_path}")
|
|
||||||
return xr.open_dataset(cache_path)
|
|
||||||
|
|
||||||
import pystac_client
|
|
||||||
import planetary_computer
|
|
||||||
import odc.stac
|
|
||||||
|
|
||||||
catalog = pystac_client.Client.open(
|
|
||||||
"https://planetarycomputer.microsoft.com/api/stac/v1",
|
|
||||||
modifier=planetary_computer.sign_inplace,
|
|
||||||
)
|
|
||||||
search = catalog.search(
|
|
||||||
collections=["sentinel-2-l2a"],
|
|
||||||
bbox=bbox,
|
|
||||||
datetime=f"{date_range[0]}/{date_range[1]}",
|
|
||||||
)
|
|
||||||
items = list(search.items())
|
|
||||||
|
|
||||||
data = odc.stac.load(
|
|
||||||
items,
|
|
||||||
bands=["red", "nir", "SCL"],
|
|
||||||
bbox=bbox,
|
|
||||||
crs="EPSG:32648",
|
|
||||||
resolution=RESOLUTION,
|
|
||||||
chunks={"x": 2048, "y": 2048, "time": 1},
|
|
||||||
groupby="solar_day"
|
|
||||||
)
|
|
||||||
if "SCL" in data.data_vars:
|
|
||||||
data = data.rename({"SCL": "scl"})
|
|
||||||
|
|
||||||
data = data.compute()
|
|
||||||
print(f"💾 Caching S2 (coord) data to {cache_path}")
|
|
||||||
data.to_netcdf(cache_path)
|
|
||||||
|
|
||||||
return data"""
|
|
||||||
content = re.sub(r'def load_data_sen2\(dc, date_range, coordinates\):.*?return data', load_data_sen2_replacement, content, flags=re.DOTALL)
|
|
||||||
|
|
||||||
|
|
||||||
with open('/home/x79/remote-sensing/new_import_ODC.py', 'w', encoding='utf-8') as f:
|
|
||||||
f.write(content)
|
|
||||||
|
|
||||||
print("Patching successful.")
|
|
||||||
@@ -1,29 +0,0 @@
|
|||||||
import json
|
|
||||||
import glob
|
|
||||||
|
|
||||||
def fix_load_sen1(file_path):
|
|
||||||
with open(file_path, 'r', encoding='utf-8') as f:
|
|
||||||
nb = json.load(f)
|
|
||||||
|
|
||||||
changed = False
|
|
||||||
for cell in nb.get('cells', []):
|
|
||||||
if cell.get('cell_type') == 'code':
|
|
||||||
source = cell.get('source', [])
|
|
||||||
for i, line in enumerate(source):
|
|
||||||
if 'load_sen1(name_vh, name_vv)' in line:
|
|
||||||
indent = line[:len(line) - len(line.lstrip())]
|
|
||||||
replacement = (
|
|
||||||
f"{indent}bbox = [longtitude_range[0], latitude_range[0], longtitude_range[1], latitude_range[1]]\n"
|
|
||||||
f"{indent}time_range = f'{{date_range[0]}}/{{date_range[1]}}'\n"
|
|
||||||
f"{indent}{line.lstrip().replace('load_sen1(name_vh, name_vv)', 'load_sen1(bbox, time_range)')}"
|
|
||||||
)
|
|
||||||
source[i] = replacement
|
|
||||||
changed = True
|
|
||||||
|
|
||||||
if changed:
|
|
||||||
with open(file_path, 'w', encoding='utf-8') as f:
|
|
||||||
json.dump(nb, f, indent=1)
|
|
||||||
print(f"Patched load_sen1 in {file_path}")
|
|
||||||
|
|
||||||
for nb in glob.glob("*.ipynb"):
|
|
||||||
fix_load_sen1(nb)
|
|
||||||
@@ -1,43 +0,0 @@
|
|||||||
import json
|
|
||||||
import glob
|
|
||||||
|
|
||||||
def patch_notebook(file_path):
|
|
||||||
with open(file_path, 'r', encoding='utf-8') as f:
|
|
||||||
nb = json.load(f)
|
|
||||||
|
|
||||||
changed = False
|
|
||||||
for cell in nb.get('cells', []):
|
|
||||||
if cell.get('cell_type') == 'code':
|
|
||||||
source = cell.get('source', [])
|
|
||||||
|
|
||||||
# Check if this cell should be fully commented out
|
|
||||||
full_source = ''.join(source)
|
|
||||||
if 'dc.load(' in full_source or 'ds.vv' in full_source:
|
|
||||||
for i in range(len(source)):
|
|
||||||
if not source[i].startswith('#'):
|
|
||||||
source[i] = '# ' + source[i]
|
|
||||||
changed = True
|
|
||||||
continue
|
|
||||||
|
|
||||||
# Otherwise, do line-by-line replacements
|
|
||||||
for i, line in enumerate(source):
|
|
||||||
if 'ST_training data_updated_1130points.shp' in line:
|
|
||||||
source[i] = line.replace('ST_training data_updated_1130points.shp', 'ST_training_data_updated_1130points.shp')
|
|
||||||
changed = True
|
|
||||||
if 'from new_import import *' in line:
|
|
||||||
source[i] = line.replace('from new_import import *', 'from new_import_ODC import *')
|
|
||||||
changed = True
|
|
||||||
if 'dc = datacube.Datacube()' in line:
|
|
||||||
source[i] = line.replace('dc = datacube.Datacube()', 'dc = None')
|
|
||||||
changed = True
|
|
||||||
if 'load_data(dc,' in line:
|
|
||||||
source[i] = line.replace('load_data(dc,', 'load_data(None,')
|
|
||||||
changed = True
|
|
||||||
|
|
||||||
if changed:
|
|
||||||
with open(file_path, 'w', encoding='utf-8') as f:
|
|
||||||
json.dump(nb, f, indent=1)
|
|
||||||
print(f"Patched {file_path}")
|
|
||||||
|
|
||||||
for nb in glob.glob("*.ipynb"):
|
|
||||||
patch_notebook(nb)
|
|
||||||
@@ -1,75 +0,0 @@
|
|||||||
import json
|
|
||||||
import glob
|
|
||||||
import re
|
|
||||||
|
|
||||||
def patch_python_script(filepath):
|
|
||||||
try:
|
|
||||||
with open(filepath, 'r', encoding='utf-8') as f:
|
|
||||||
content = f.read()
|
|
||||||
|
|
||||||
original_content = content
|
|
||||||
|
|
||||||
# Replace imports
|
|
||||||
content = re.sub(r'from sklearn\.ensemble import RandomForestClassifier',
|
|
||||||
'from xgboost import XGBClassifier', content)
|
|
||||||
|
|
||||||
# Replace the model instantiations (for 01.train_ODC.py)
|
|
||||||
rf_pattern = re.compile(r'model\s*=\s*RandomForestClassifier\([^)]+\)', re.DOTALL)
|
|
||||||
xgb_replacement = """model = XGBClassifier(
|
|
||||||
n_estimators=200,
|
|
||||||
max_depth=30,
|
|
||||||
tree_method="hist",
|
|
||||||
device="cuda",
|
|
||||||
random_state=42,
|
|
||||||
n_jobs=-1,
|
|
||||||
verbosity=1
|
|
||||||
)"""
|
|
||||||
|
|
||||||
content = rf_pattern.sub(xgb_replacement, content)
|
|
||||||
|
|
||||||
if content != original_content:
|
|
||||||
with open(filepath, 'w', encoding='utf-8') as f:
|
|
||||||
f.write(content)
|
|
||||||
print(f"Patched {filepath}")
|
|
||||||
except Exception as e:
|
|
||||||
print(f"Error patching {filepath}: {e}")
|
|
||||||
|
|
||||||
def patch_notebook(filepath):
|
|
||||||
try:
|
|
||||||
with open(filepath, 'r', encoding='utf-8') as f:
|
|
||||||
nb = json.load(f)
|
|
||||||
|
|
||||||
changed = False
|
|
||||||
import_pattern = re.compile(r'from\s+sklearn\.ensemble\s+import\s+RandomForestClassifier')
|
|
||||||
inst_pattern = re.compile(r'RandomForestClassifier\([^)]*\)')
|
|
||||||
|
|
||||||
for cell in nb.get('cells', []):
|
|
||||||
if cell.get('cell_type') == 'code':
|
|
||||||
source = cell.get('source', [])
|
|
||||||
for i in range(len(source)):
|
|
||||||
if import_pattern.search(source[i]):
|
|
||||||
source[i] = import_pattern.sub('from xgboost import XGBClassifier', source[i])
|
|
||||||
changed = True
|
|
||||||
|
|
||||||
if inst_pattern.search(source[i]):
|
|
||||||
source[i] = inst_pattern.sub("XGBClassifier(tree_method='hist', device='cuda', random_state=42, n_jobs=-1)", source[i])
|
|
||||||
changed = True
|
|
||||||
|
|
||||||
if "'classifier__criterion': ['gini', 'entropy']" in source[i]:
|
|
||||||
source[i] = source[i].replace("'classifier__criterion': ['gini', 'entropy']",
|
|
||||||
"'classifier__learning_rate': [0.01, 0.1, 0.2]")
|
|
||||||
changed = True
|
|
||||||
|
|
||||||
if changed:
|
|
||||||
with open(filepath, 'w', encoding='utf-8') as f:
|
|
||||||
json.dump(nb, f, indent=1)
|
|
||||||
print(f"Patched {filepath}")
|
|
||||||
except Exception as e:
|
|
||||||
print(f"Error patching {filepath}: {e}")
|
|
||||||
|
|
||||||
if __name__ == "__main__":
|
|
||||||
patch_python_script("01.train_ODC.py")
|
|
||||||
patch_python_script("new_train.py")
|
|
||||||
|
|
||||||
for nb in glob.glob("*.ipynb"):
|
|
||||||
patch_notebook(nb)
|
|
||||||
@@ -1,272 +0,0 @@
|
|||||||
#!/usr/bin/env python
|
|
||||||
# coding: utf-8
|
|
||||||
|
|
||||||
# In[1]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '%matplotlib inline\nfrom new_import import *\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[2]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Dask gateway\ncluster, client = notebook_utils.initialize_dask(use_gateway=True, workers=(1,4))\ndc = datacube.Datacube()\n\n# Configure s3 access\nconfigure_s3_access(aws_unsigned=False, requester_pays=True, client=client)\n\nclient\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[3]:
|
|
||||||
|
|
||||||
|
|
||||||
## cấu hình thời gian lấy ảnh và tọa độ
|
|
||||||
# date_range = ('2022-09-01', '2023-10-01')
|
|
||||||
# longtitude_range = (105.86575, 105.94120)
|
|
||||||
# latitude_range = (9.65070, 9.69850)
|
|
||||||
date_range = ('2022-09-01', '2023-10-01')
|
|
||||||
longtitude_range = (105.5, 106.4)
|
|
||||||
latitude_range = (9.2, 10.0)
|
|
||||||
|
|
||||||
|
|
||||||
# In[4]:
|
|
||||||
|
|
||||||
|
|
||||||
## truy vấn ảnh vệ tinh sen2
|
|
||||||
data = load_data(dc, date_range, longtitude_range, latitude_range)
|
|
||||||
notebook_utils.heading(notebook_utils.xarray_object_size(data))
|
|
||||||
display(data)
|
|
||||||
|
|
||||||
|
|
||||||
# In[5]:
|
|
||||||
|
|
||||||
|
|
||||||
# Specify the start and end times
|
|
||||||
min_date = '2022-09-01' # Thời gian bắt đầu lấy data cho quá trình train
|
|
||||||
max_date = '2023-10-01' # Thời gian kết thúc lấy data cho quá trình train
|
|
||||||
# Just do 1 month for testing
|
|
||||||
# max_date = '2022-10-01' # Thời gian kết thúc lấy data cho quá trình train
|
|
||||||
|
|
||||||
# Specify a spatail region to search using latitude/longitude cooridinates
|
|
||||||
min_longitude, max_longitude = (105.5, 106.4)
|
|
||||||
min_latitude, max_latitude = (9.2, 10.0)
|
|
||||||
|
|
||||||
# Specify the product. In this case we want to use Sentinel-2 Level-2A data
|
|
||||||
product = 's2_l2a'
|
|
||||||
|
|
||||||
# Construct the search query dictionary
|
|
||||||
query = {
|
|
||||||
'product': product, # Product name
|
|
||||||
'x': (min_longitude, max_longitude), # "x" axis bounds
|
|
||||||
'y': (min_latitude, max_latitude), # "y" axis bounds
|
|
||||||
'time': (min_date, max_date), # Any parsable date strings
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
# In[6]:
|
|
||||||
|
|
||||||
|
|
||||||
# Most common CRS
|
|
||||||
native_crs = notebook_utils.mostcommon_crs(dc, query)
|
|
||||||
print(f'Most common native CRS: {native_crs}')
|
|
||||||
|
|
||||||
|
|
||||||
# In[7]:
|
|
||||||
|
|
||||||
|
|
||||||
# Specify the spectral band measurements we want to use for a classification algorithm
|
|
||||||
measurements = ['red', 'nir', 'scl']
|
|
||||||
|
|
||||||
load_params = {
|
|
||||||
'measurements': measurements, # Selected measurement or alias names
|
|
||||||
'output_crs': native_crs, # Target EPSG code
|
|
||||||
'resolution': (-10, 10), # Target resolution
|
|
||||||
'group_by': 'solar_day', # Scene grouping
|
|
||||||
'dask_chunks': {'x': 2048, 'y': 2048}, # Dask chunks
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
# In[8]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# The replacement "dc.load()" function for this product\ndata = load_s2l2a_with_offset(\n dc,\n query | load_params # Combine the two dicts that contain our search and load parameters\n)\n\n# This line prints the total size of the dataset hat was loaded\nnotebook_utils.heading(notebook_utils.xarray_object_size(data))\n\ndisplay(data)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[9]:
|
|
||||||
|
|
||||||
|
|
||||||
# %%time
|
|
||||||
# # Tiến hành loại bỏ các vị trí bị mây ảnh hưởng
|
|
||||||
# result = mask_clean(data)
|
|
||||||
# progress(result)
|
|
||||||
|
|
||||||
|
|
||||||
# In[10]:
|
|
||||||
|
|
||||||
|
|
||||||
# Tiến hành tính toán NDVI
|
|
||||||
ds1 = calculate_indices(result, index='NDVI', satellite_mission='s2')
|
|
||||||
ndvi = ds1["NDVI"]
|
|
||||||
display(ndvi)
|
|
||||||
|
|
||||||
|
|
||||||
# In[11]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', "## tính ndvi theo tháng\naverage_ndvi = ndvi.resample(time='1M').mean().persist()\nprogress(average_ndvi)\n")
|
|
||||||
|
|
||||||
|
|
||||||
# In[12]:
|
|
||||||
|
|
||||||
|
|
||||||
# compute average_ndvi
|
|
||||||
average_ndvi = average_ndvi.compute()
|
|
||||||
|
|
||||||
|
|
||||||
# In[13]:
|
|
||||||
|
|
||||||
|
|
||||||
# cấu hình vh vv file
|
|
||||||
# name_vh = "ThuanHoa/ThuanHoa_VH.tif"
|
|
||||||
# name_vv = "ThuanHoa/ThuanHoa_VV.tif"
|
|
||||||
|
|
||||||
# load dữ liệu sen1
|
|
||||||
bbox = [105.5, 9.2, 106.4, 10.0]
|
|
||||||
time_range = '2022-09-01/2023-10-01'
|
|
||||||
# dsvh, dsvv = load_sen1(bbox, time_range)
|
|
||||||
|
|
||||||
name_vh = "vh-0922_0923-full_ST.tif"
|
|
||||||
name_vv = "vv-0922_0923-full_ST.tif"
|
|
||||||
|
|
||||||
if not os.path.exists(name_vh):
|
|
||||||
get_ipython().system('aws s3 cp s3://easi-asia-dc-data/staging/ctu/sentinel-1/vh-0922_0923-full_ST.tif vh-0922_0923-full_ST.tif')
|
|
||||||
if not os.path.exists(name_vv):
|
|
||||||
get_ipython().system('aws s3 cp s3://easi-asia-dc-data/staging/ctu/sentinel-1/vv-0922_0923-full_ST.tif vv-0922_0923-full_ST.tif')
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.4, 10.0]
|
|
||||||
time_range = '2022-09-01/2023-10-01'
|
|
||||||
dsvh, dsvv = load_sen1(bbox, time_range)
|
|
||||||
|
|
||||||
|
|
||||||
# In[27]:
|
|
||||||
|
|
||||||
|
|
||||||
from sklearn.preprocessing import PolynomialFeatures
|
|
||||||
from sklearn.linear_model import LinearRegression
|
|
||||||
from sklearn.ensemble import RandomForestRegressor
|
|
||||||
|
|
||||||
|
|
||||||
# In[28]:
|
|
||||||
|
|
||||||
|
|
||||||
average_ndvi = average_ndvi[:, :7680, :8687]
|
|
||||||
mask = ~np.isnan(average_ndvi)
|
|
||||||
print(average_ndvi.shape)
|
|
||||||
print(dsvh.shape)
|
|
||||||
print(dsvv.shape)
|
|
||||||
print(mask.shape)
|
|
||||||
X_train = np.stack([dsvh.values[mask], dsvv.values[mask]], axis=1)
|
|
||||||
y_train = average_ndvi.values[mask]
|
|
||||||
|
|
||||||
|
|
||||||
# In[29]:
|
|
||||||
|
|
||||||
|
|
||||||
model = LinearRegression()
|
|
||||||
model.fit(X_train, y_train)
|
|
||||||
|
|
||||||
|
|
||||||
# In[30]:
|
|
||||||
|
|
||||||
|
|
||||||
X_pred = np.stack([dsvh.values[~mask], dsvv.values[~mask]], axis=1)
|
|
||||||
average_ndvi.values[~mask] = model.predict(X_pred)
|
|
||||||
|
|
||||||
|
|
||||||
# In[31]:
|
|
||||||
|
|
||||||
|
|
||||||
average_ndvi_filled = xr.DataArray(average_ndvi, dims=average_ndvi.dims)
|
|
||||||
|
|
||||||
|
|
||||||
# In[32]:
|
|
||||||
|
|
||||||
|
|
||||||
plt.imshow(average_ndvi_filled.isel(time=6))
|
|
||||||
|
|
||||||
|
|
||||||
# In[65]:
|
|
||||||
|
|
||||||
|
|
||||||
plt.imshow(average_ndvi.isel(time=6))
|
|
||||||
|
|
||||||
|
|
||||||
# In[33]:
|
|
||||||
|
|
||||||
|
|
||||||
train_path = "train/ST_training data_updated_1130points.shp"
|
|
||||||
|
|
||||||
|
|
||||||
# In[34]:
|
|
||||||
|
|
||||||
|
|
||||||
train = load_train_data(train_path)
|
|
||||||
|
|
||||||
|
|
||||||
# In[37]:
|
|
||||||
|
|
||||||
|
|
||||||
datasets = get_data_sen1_and_sen2(train, average_ndvi_filled, dsvh, dsvv)
|
|
||||||
|
|
||||||
|
|
||||||
# In[39]:
|
|
||||||
|
|
||||||
|
|
||||||
# cấu hình nhãn dữ liệu
|
|
||||||
label_mapping = {
|
|
||||||
"Lua tom": "0",
|
|
||||||
"Lua": "1",
|
|
||||||
"CHN": "2",
|
|
||||||
"CLN": "3",
|
|
||||||
"TS": "4",
|
|
||||||
"Song": "5",
|
|
||||||
"Dat xay dung": "6",
|
|
||||||
"Rung": "7"
|
|
||||||
}
|
|
||||||
|
|
||||||
# chia tập dữ liệu train, val, test
|
|
||||||
X_train, X_val, X_test, y_train, y_val, y_test = split_train_data(train, label_mapping, datasets)
|
|
||||||
|
|
||||||
|
|
||||||
# In[40]:
|
|
||||||
|
|
||||||
|
|
||||||
# Huấn luyện mô hình
|
|
||||||
grid_search = train_with_rf(X_train, X_val, y_train, y_val)
|
|
||||||
|
|
||||||
|
|
||||||
# In[41]:
|
|
||||||
|
|
||||||
|
|
||||||
# kiểm tra độ chính xác với tập test
|
|
||||||
y_pred_test = grid_search.predict(X_test)
|
|
||||||
test_accuracy = accuracy_score(y_test, y_pred_test)
|
|
||||||
print(f"Accuracy for test data {round(test_accuracy, 2)*100} %")
|
|
||||||
|
|
||||||
|
|
||||||
# In[42]:
|
|
||||||
|
|
||||||
|
|
||||||
# Lưu mô hình huấn luyện
|
|
||||||
save_model("model_new.joblib", grid_search)
|
|
||||||
|
|
||||||
|
|
||||||
# In[43]:
|
|
||||||
|
|
||||||
|
|
||||||
# đóng client, cluster
|
|
||||||
client.close()
|
|
||||||
cluster.close()
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
@@ -1,127 +0,0 @@
|
|||||||
#!/usr/bin/env python
|
|
||||||
# coding: utf-8
|
|
||||||
|
|
||||||
# In[1]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '%matplotlib inline\nfrom new_import import *\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[2]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Dask gateway\ncluster, client = notebook_utils.initialize_dask(use_gateway=True, workers=(1,4))\ndc = datacube.Datacube()\n\n# Configure s3 access\nconfigure_s3_access(aws_unsigned=False, requester_pays=True, client=client)\n\nclient\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[3]:
|
|
||||||
|
|
||||||
|
|
||||||
## cấu hình thời gian lấy ảnh và tọa độ
|
|
||||||
date_range = ('2022-09-01', '2023-10-01')
|
|
||||||
longtitude_range = (105.86575, 105.94120)
|
|
||||||
latitude_range = (9.65070, 9.69850)
|
|
||||||
|
|
||||||
|
|
||||||
# In[4]:
|
|
||||||
|
|
||||||
|
|
||||||
## truy vấn ảnh vệ tinh sen2
|
|
||||||
data = load_data(dc, date_range, longtitude_range, latitude_range)
|
|
||||||
notebook_utils.heading(notebook_utils.xarray_object_size(data))
|
|
||||||
display(data)
|
|
||||||
|
|
||||||
|
|
||||||
# In[5]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', '# Tiến hành loại bỏ các vị trí bị mây ảnh hưởng\nresult = mask_clean(data)\nprogress(result)\n')
|
|
||||||
|
|
||||||
|
|
||||||
# In[6]:
|
|
||||||
|
|
||||||
|
|
||||||
# Tiến hành tính toán NDVI
|
|
||||||
ds1 = calculate_indices(result, index='NDVI', satellite_mission='s2')
|
|
||||||
ndvi = ds1["NDVI"]
|
|
||||||
display(ndvi)
|
|
||||||
|
|
||||||
|
|
||||||
# In[17]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_cell_magic('time', '', "## tính ndvi theo tháng\naverage_ndvi = ndvi.resample(time='1M').mean().persist()\nprogress(average_ndvi)\n")
|
|
||||||
|
|
||||||
|
|
||||||
# In[18]:
|
|
||||||
|
|
||||||
|
|
||||||
# compute average_ndvi
|
|
||||||
average_ndvi = average_ndvi.compute()
|
|
||||||
|
|
||||||
|
|
||||||
# In[9]:
|
|
||||||
|
|
||||||
|
|
||||||
# cấu hình vh vv file
|
|
||||||
name_vh = "ThuanHoa/ThuanHoa_VH.tif"
|
|
||||||
name_vv = "ThuanHoa/ThuanHoa_VV.tif"
|
|
||||||
|
|
||||||
# load dữ liệu sen1
|
|
||||||
bbox = [105.5, 9.2, 106.4, 10.0]
|
|
||||||
time_range = '2022-09-01/2023-10-01'
|
|
||||||
dsvh, dsvv = load_sen1(bbox, time_range)
|
|
||||||
|
|
||||||
|
|
||||||
# In[10]:
|
|
||||||
|
|
||||||
|
|
||||||
from sklearn.preprocessing import PolynomialFeatures
|
|
||||||
from sklearn.linear_model import LinearRegression
|
|
||||||
from sklearn.ensemble import RandomForestRegressor
|
|
||||||
|
|
||||||
|
|
||||||
# In[11]:
|
|
||||||
|
|
||||||
|
|
||||||
mask = ~np.isnan(average_ndvi)
|
|
||||||
X_train = np.stack([dsvh.values[mask], dsvv.values[mask]], axis=1)
|
|
||||||
y_train = average_ndvi.values[mask]
|
|
||||||
|
|
||||||
|
|
||||||
# In[12]:
|
|
||||||
|
|
||||||
|
|
||||||
model = LinearRegression()
|
|
||||||
model.fit(X_train, y_train)
|
|
||||||
|
|
||||||
|
|
||||||
# In[13]:
|
|
||||||
|
|
||||||
|
|
||||||
X_pred = np.stack([dsvh.values[~mask], dsvv.values[~mask]], axis=1)
|
|
||||||
average_ndvi.values[~mask] = model.predict(X_pred)
|
|
||||||
|
|
||||||
|
|
||||||
# In[14]:
|
|
||||||
|
|
||||||
|
|
||||||
average_ndvi_filled = xr.DataArray(average_ndvi, dims=average_ndvi.dims)
|
|
||||||
|
|
||||||
|
|
||||||
# In[16]:
|
|
||||||
|
|
||||||
|
|
||||||
plt.imshow(average_ndvi_filled.isel(time=1))
|
|
||||||
|
|
||||||
|
|
||||||
# In[19]:
|
|
||||||
|
|
||||||
|
|
||||||
plt.imshow(average_ndvi.isel(time=1))
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
@@ -1,145 +0,0 @@
|
|||||||
import json
|
|
||||||
|
|
||||||
def get_content(filename):
|
|
||||||
with open(filename, "r", encoding="utf-8") as f:
|
|
||||||
return f.read()
|
|
||||||
|
|
||||||
fe_content = get_content("feature_extractor.py")
|
|
||||||
tm_content = get_content("train_module.py")
|
|
||||||
dt_content = get_content("train_land_decision_tree_gpu.py")
|
|
||||||
|
|
||||||
notebook = {
|
|
||||||
"cells": [
|
|
||||||
{
|
|
||||||
"cell_type": "markdown",
|
|
||||||
"metadata": {},
|
|
||||||
"source": [
|
|
||||||
"# Tải Dữ liệu Vệ tinh qua Colab (Self-contained)\n",
|
|
||||||
"Notebook này đã được nhúng sẵn toàn bộ mã nguồn xử lý. Bạn không cần upload cả thư mục `remote-sensing` nữa.\n",
|
|
||||||
"\n",
|
|
||||||
"## Bước 1: Upload Shapefile (BẮT BUỘC)\n",
|
|
||||||
"Mô hình cần biết các điểm tọa độ đất để lấy dữ liệu. Hãy nén thư mục `train/` trên máy bạn thành `train.zip` và chạy ô dưới đây để upload nó trực tiếp lên Colab."
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "code",
|
|
||||||
"execution_count": None,
|
|
||||||
"metadata": {},
|
|
||||||
"outputs": [],
|
|
||||||
"source": [
|
|
||||||
"from google.colab import files\n",
|
|
||||||
"import os\n",
|
|
||||||
"\n",
|
|
||||||
"print(\"Hãy chọn file train.zip từ máy tính của bạn:\")\n",
|
|
||||||
"uploaded = files.upload()\n",
|
|
||||||
"\n",
|
|
||||||
"if \"train.zip\" in uploaded:\n",
|
|
||||||
" !unzip -q -o train.zip -d /content/train_tmp/\n",
|
|
||||||
" # Move the extracted files directly to /content/train/\n",
|
|
||||||
" !mkdir -p /content/train\n",
|
|
||||||
" !mv /content/train_tmp/*/* /content/train/ 2>/dev/null || mv /content/train_tmp/* /content/train/\n",
|
|
||||||
" print(\"Đã giải nén shapefile thành công vào thư mục /content/train/\")\n",
|
|
||||||
"else:\n",
|
|
||||||
" print(\"LỖI: Bạn chưa upload file có tên là train.zip!\")"
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "markdown",
|
|
||||||
"metadata": {},
|
|
||||||
"source": [
|
|
||||||
"## Bước 2: Cài đặt thư viện & Tạo môi trường"
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "code",
|
|
||||||
"execution_count": None,
|
|
||||||
"metadata": {},
|
|
||||||
"outputs": [],
|
|
||||||
"source": [
|
|
||||||
"!pip install planetary-computer pystac-client odc-stac geopandas rasterio xarray joblib scikit-learn xgboost lightgbm"
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "code",
|
|
||||||
"execution_count": None,
|
|
||||||
"metadata": {},
|
|
||||||
"outputs": [],
|
|
||||||
"source": [
|
|
||||||
"%%writefile feature_extractor.py\n" + fe_content
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "code",
|
|
||||||
"execution_count": None,
|
|
||||||
"metadata": {},
|
|
||||||
"outputs": [],
|
|
||||||
"source": [
|
|
||||||
"%%writefile train_module.py\n" + tm_content
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "code",
|
|
||||||
"execution_count": None,
|
|
||||||
"metadata": {},
|
|
||||||
"outputs": [],
|
|
||||||
"source": [
|
|
||||||
"%%writefile train_land_decision_tree_gpu.py\n" + dt_content
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "markdown",
|
|
||||||
"metadata": {},
|
|
||||||
"source": [
|
|
||||||
"## Bước 3: Chạy tiến trình tải ảnh vệ tinh và tạo Cache"
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "code",
|
|
||||||
"execution_count": None,
|
|
||||||
"metadata": {},
|
|
||||||
"outputs": [],
|
|
||||||
"source": [
|
|
||||||
"!python train_land_decision_tree_gpu.py"
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "markdown",
|
|
||||||
"metadata": {},
|
|
||||||
"source": [
|
|
||||||
"## Bước 4: Tải file Cache về máy"
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "code",
|
|
||||||
"execution_count": None,
|
|
||||||
"metadata": {},
|
|
||||||
"outputs": [],
|
|
||||||
"source": [
|
|
||||||
"from google.colab import files\n",
|
|
||||||
"import glob\n",
|
|
||||||
"\n",
|
|
||||||
"cache_files = glob.glob(\"dataset_cache/*.joblib\")\n",
|
|
||||||
"if cache_files:\n",
|
|
||||||
" latest_cache = max(cache_files, key=os.path.getctime)\n",
|
|
||||||
" print(f\"Đang tải file {latest_cache} về máy...\")\n",
|
|
||||||
" files.download(latest_cache)\n",
|
|
||||||
"else:\n",
|
|
||||||
" print(\"Chưa tìm thấy file cache. Hãy chắc chắn bước 3 đã chạy thành công!\")"
|
|
||||||
]
|
|
||||||
}
|
|
||||||
],
|
|
||||||
"metadata": {
|
|
||||||
"kernelspec": {
|
|
||||||
"display_name": "Python 3",
|
|
||||||
"language": "python",
|
|
||||||
"name": "python3"
|
|
||||||
}
|
|
||||||
},
|
|
||||||
"nbformat": 4,
|
|
||||||
"nbformat_minor": 4
|
|
||||||
}
|
|
||||||
|
|
||||||
with open("Download_Cache_Colab.ipynb", "w", encoding="utf-8") as f:
|
|
||||||
json.dump(notebook, f, indent=1, ensure_ascii=False)
|
|
||||||
|
|
||||||
print("Notebook updated successfully!")
|
|
||||||
Binary file not shown.
@@ -1,35 +0,0 @@
|
|||||||
import json
|
|
||||||
import glob
|
|
||||||
|
|
||||||
def fix_notebook(file_path):
|
|
||||||
with open(file_path, 'r', encoding='utf-8') as f:
|
|
||||||
nb = json.load(f)
|
|
||||||
changed = False
|
|
||||||
for cell in nb.get('cells', []):
|
|
||||||
if cell.get('cell_type') == 'code':
|
|
||||||
source = cell.get('source', [])
|
|
||||||
if isinstance(source, list):
|
|
||||||
for i, line in enumerate(source):
|
|
||||||
if "dc = datacube.Datacube()" in line:
|
|
||||||
source[i] = "dc = None\n"
|
|
||||||
changed = True
|
|
||||||
if "ds = dc.load(" in line:
|
|
||||||
source[i] = "ds = None\n"
|
|
||||||
changed = True
|
|
||||||
if "data = dc.load(" in line:
|
|
||||||
source[i] = "data = None\n"
|
|
||||||
changed = True
|
|
||||||
# If ds is None, ds.vv will fail
|
|
||||||
if "vv_data = ds.vv" in line:
|
|
||||||
source[i] = "vv_data = None\n"
|
|
||||||
changed = True
|
|
||||||
if "notebook_utils.xarray_object_size(ds)" in line:
|
|
||||||
source[i] = line.replace("notebook_utils.xarray_object_size(ds)", "'ds is None'")
|
|
||||||
changed = True
|
|
||||||
if changed:
|
|
||||||
with open(file_path, 'w', encoding='utf-8') as f:
|
|
||||||
json.dump(nb, f, indent=1)
|
|
||||||
print(f"Removed datacube from {file_path}")
|
|
||||||
|
|
||||||
for nb in glob.glob("*.ipynb"):
|
|
||||||
fix_notebook(nb)
|
|
||||||
@@ -1,110 +0,0 @@
|
|||||||
import os
|
|
||||||
import glob
|
|
||||||
import json
|
|
||||||
import subprocess
|
|
||||||
import time
|
|
||||||
from tabulate import tabulate
|
|
||||||
|
|
||||||
scripts = [
|
|
||||||
"train_land_randomforest.py",
|
|
||||||
"train_cloud_cnn.py",
|
|
||||||
"train_cloud_swin_unet.py",
|
|
||||||
"train_ndvi_statistical.py",
|
|
||||||
"train_ndvi_lstm_gru.py",
|
|
||||||
"train_ndvi_convlstm.py",
|
|
||||||
"train_ndvi_hybrid_physics.py",
|
|
||||||
"train_ndvi_ensemble.py"
|
|
||||||
]
|
|
||||||
|
|
||||||
print("🚀 Đang khởi chạy song song tất cả các mô hình...")
|
|
||||||
processes = []
|
|
||||||
for script in scripts:
|
|
||||||
if os.path.exists(script):
|
|
||||||
cmd = f"source /home/x79/miniconda3/etc/profile.d/conda.sh && conda activate env_01 && python {script}"
|
|
||||||
p = subprocess.Popen(["bash", "-c", cmd], stdout=subprocess.DEVNULL, stderr=subprocess.DEVNULL)
|
|
||||||
processes.append((script, p))
|
|
||||||
|
|
||||||
for script, p in processes:
|
|
||||||
p.wait()
|
|
||||||
|
|
||||||
print("✅ Đã chạy xong tất cả các mô hình!\n")
|
|
||||||
print("📊 BẢNG SO SÁNH KẾT QUẢ CÁC MÔ HÌNH\n")
|
|
||||||
|
|
||||||
# 1. Phân loại đất
|
|
||||||
print("### 1. Nhóm Phân loại Lớp phủ (Land Classification)")
|
|
||||||
land_data = []
|
|
||||||
|
|
||||||
# Đọc XGBoost từ thư mục gốc
|
|
||||||
if os.path.exists("model_xgboost_info.json"):
|
|
||||||
with open("model_xgboost_info.json", 'r') as f:
|
|
||||||
data = json.load(f)
|
|
||||||
params = data.get('params', {})
|
|
||||||
param_str = f"estimators:{params.get('n_estimators')}, depth:{params.get('max_depth')}" if params else "N/A"
|
|
||||||
land_data.append([
|
|
||||||
data.get('model_type', 'XGBoost'),
|
|
||||||
data.get('accuracy', ''),
|
|
||||||
data.get('precision', ''),
|
|
||||||
data.get('recall', ''),
|
|
||||||
data.get('f1_score', ''),
|
|
||||||
param_str
|
|
||||||
])
|
|
||||||
|
|
||||||
# Đọc các model khác trong model_train
|
|
||||||
for info_file in glob.glob("model_train/*_info.json"):
|
|
||||||
with open(info_file, 'r') as f:
|
|
||||||
data = json.load(f)
|
|
||||||
# Chỉ lấy các model có độ chính xác (để lọc model rác/cũ)
|
|
||||||
if 'accuracy' not in data and 'f1_score' not in data:
|
|
||||||
continue
|
|
||||||
|
|
||||||
params = data.get('params', {})
|
|
||||||
param_str = f"estimators:{params.get('n_estimators')}, depth:{params.get('max_depth')}" if params else "N/A"
|
|
||||||
|
|
||||||
# Fallback for Random Forest
|
|
||||||
if data.get('model_type') == 'RandomForest_RealData':
|
|
||||||
param_str = "estimators:100, depth:15"
|
|
||||||
|
|
||||||
land_data.append([
|
|
||||||
data.get('model_type', ''),
|
|
||||||
data.get('accuracy', ''),
|
|
||||||
data.get('precision', ''),
|
|
||||||
data.get('recall', ''),
|
|
||||||
data.get('f1_score', ''),
|
|
||||||
param_str
|
|
||||||
])
|
|
||||||
|
|
||||||
if land_data:
|
|
||||||
print(tabulate(land_data, headers=["Model", "Accuracy", "Precision", "Recall", "F1-Score", "Parameters"], tablefmt="github"))
|
|
||||||
print("\n")
|
|
||||||
|
|
||||||
# 2. Xóa mây
|
|
||||||
print("### 2. Nhóm Xóa mây (Cloud Removal)")
|
|
||||||
cloud_data = []
|
|
||||||
for info_file in glob.glob("cloud_removal_model/*_info.json"):
|
|
||||||
with open(info_file, 'r') as f:
|
|
||||||
data = json.load(f)
|
|
||||||
cloud_data.append([
|
|
||||||
data.get('model_type', ''),
|
|
||||||
data.get('epoch', ''),
|
|
||||||
data.get('train_loss', ''),
|
|
||||||
data.get('val_loss', '')
|
|
||||||
])
|
|
||||||
if cloud_data:
|
|
||||||
print(tabulate(cloud_data, headers=["Model", "Epochs", "Train Loss", "Val Loss"], tablefmt="github"))
|
|
||||||
print("\n")
|
|
||||||
|
|
||||||
# 3. Dự báo NDVI
|
|
||||||
print("### 3. Nhóm Dự báo Thực vật (NDVI Forecasting)")
|
|
||||||
ndvi_data = []
|
|
||||||
for info_file in glob.glob("ndvi_forecast_model/*_info.json"):
|
|
||||||
with open(info_file, 'r') as f:
|
|
||||||
data = json.load(f)
|
|
||||||
ndvi_data.append([
|
|
||||||
data.get('model_type', ''),
|
|
||||||
data.get('rmse', ''),
|
|
||||||
data.get('mae', ''),
|
|
||||||
data.get('epoch', 'N/A')
|
|
||||||
])
|
|
||||||
if ndvi_data:
|
|
||||||
print(tabulate(ndvi_data, headers=["Model", "RMSE", "MAE", "Epochs"], tablefmt="github"))
|
|
||||||
print("\n")
|
|
||||||
@@ -1,72 +0,0 @@
|
|||||||
#!/usr/bin/env python
|
|
||||||
# coding: utf-8
|
|
||||||
|
|
||||||
# In[2]:
|
|
||||||
|
|
||||||
|
|
||||||
get_ipython().run_line_magic('matplotlib', 'inline')
|
|
||||||
from new_import import *
|
|
||||||
|
|
||||||
|
|
||||||
# Dask gateway
|
|
||||||
cluster, client = notebook_utils.initialize_dask(use_gateway=True, workers=(1,4))
|
|
||||||
dc = datacube.Datacube()
|
|
||||||
|
|
||||||
|
|
||||||
# Configure s3 access
|
|
||||||
configure_s3_access(aws_unsigned=False, requester_pays=True, client=client)
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
# In[3]:
|
|
||||||
|
|
||||||
|
|
||||||
ds = dc.load(
|
|
||||||
product="sentinel1_grd_gamma0_20m",
|
|
||||||
x=(105.5, 106.4),
|
|
||||||
y=(9.2, 10.0),
|
|
||||||
time=("2022-09-01", "2023-10-01"),
|
|
||||||
measurements=["vv", "vh"],
|
|
||||||
output_crs="EPSG:32648",
|
|
||||||
resolution=(-10,10),
|
|
||||||
dask_chunks={"x":2048, "y":2048},
|
|
||||||
skip_broken_datasets=True,
|
|
||||||
group_by="solar_day"
|
|
||||||
)
|
|
||||||
notebook_utils.heading(notebook_utils.xarray_object_size(ds))
|
|
||||||
ds
|
|
||||||
|
|
||||||
|
|
||||||
# In[18]:
|
|
||||||
|
|
||||||
|
|
||||||
vh = ds.vh.resample(time='1M').mean().persist()
|
|
||||||
vh = vh.compute()
|
|
||||||
vv = ds.vv.resample(time='1M').mean().persist()
|
|
||||||
vv = vv.compute()
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
# In[28]:
|
|
||||||
|
|
||||||
|
|
||||||
vv.min()
|
|
||||||
|
|
||||||
|
|
||||||
# In[33]:
|
|
||||||
|
|
||||||
|
|
||||||
import matplotlib.pyplot as plt
|
|
||||||
|
|
||||||
# Plot the data
|
|
||||||
plt.imshow(vh.isel(time=0), cmap='viridis', vmin=0, vmax=1)
|
|
||||||
plt.colorbar() # Add colorbar for reference
|
|
||||||
plt.show()
|
|
||||||
|
|
||||||
|
|
||||||
# In[ ]:
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
@@ -1,33 +0,0 @@
|
|||||||
import torch
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
device = "cpu"
|
|
||||||
input_array = np.zeros((4, 16, 16), dtype=np.float32)
|
|
||||||
input_tensor = torch.from_numpy(input_array).unsqueeze(0).to(device)
|
|
||||||
|
|
||||||
print("Before pad:", input_tensor.shape)
|
|
||||||
|
|
||||||
from train_cloud_removal import UNet
|
|
||||||
model = UNet(in_channels=6, out_channels=4).to(device)
|
|
||||||
|
|
||||||
if hasattr(model, 'inc') and hasattr(model.inc.double_conv[0], 'in_channels'):
|
|
||||||
expected_channels = model.inc.double_conv[0].in_channels
|
|
||||||
elif hasattr(model, 'conv1') and hasattr(model.conv1, 'in_channels'):
|
|
||||||
expected_channels = model.conv1.in_channels
|
|
||||||
else:
|
|
||||||
expected_channels = list(model.parameters())[0].shape[1]
|
|
||||||
|
|
||||||
print("Expected channels:", expected_channels)
|
|
||||||
|
|
||||||
if expected_channels > input_tensor.shape[1]:
|
|
||||||
pad_channels = expected_channels - input_tensor.shape[1]
|
|
||||||
padding = torch.zeros(1, pad_channels, *input_tensor.shape[2:]).to(device)
|
|
||||||
input_tensor = torch.cat([input_tensor, padding], dim=1)
|
|
||||||
|
|
||||||
print("After pad:", input_tensor.shape)
|
|
||||||
|
|
||||||
try:
|
|
||||||
model(input_tensor)
|
|
||||||
print("Success!")
|
|
||||||
except Exception as e:
|
|
||||||
print("Error:", e)
|
|
||||||
@@ -1,25 +0,0 @@
|
|||||||
import os
|
|
||||||
import sys
|
|
||||||
sys.path.insert(0, os.getcwd())
|
|
||||||
import new_import_ODC
|
|
||||||
import time
|
|
||||||
import xarray as xr
|
|
||||||
|
|
||||||
# Mock minimal params to test load_data
|
|
||||||
date_range = ('2023-01-01', '2023-01-31')
|
|
||||||
longtitude_range = (105.0, 105.1)
|
|
||||||
latitude_range = (9.5, 9.6)
|
|
||||||
|
|
||||||
print("--- First Call (Downloading & Caching) ---")
|
|
||||||
start = time.time()
|
|
||||||
data1 = new_import_ODC.load_data(None, date_range, longtitude_range, latitude_range)
|
|
||||||
end = time.time()
|
|
||||||
print(f"Time taken: {end - start:.2f}s")
|
|
||||||
|
|
||||||
print("--- Second Call (Loading from Cache) ---")
|
|
||||||
start = time.time()
|
|
||||||
data2 = new_import_ODC.load_data(None, date_range, longtitude_range, latitude_range)
|
|
||||||
end = time.time()
|
|
||||||
print(f"Time taken: {end - start:.2f}s")
|
|
||||||
|
|
||||||
print("✅ Test completed")
|
|
||||||
@@ -1,7 +0,0 @@
|
|||||||
import torch
|
|
||||||
checkpoint = torch.load('cloud_removal_model/cloud_removal_unet_best.pth', map_location='cpu')
|
|
||||||
print(checkpoint.keys())
|
|
||||||
print("in_channels in checkpoint:", 'in_channels' in checkpoint)
|
|
||||||
if 'in_channels' in checkpoint:
|
|
||||||
print(checkpoint['in_channels'])
|
|
||||||
print("Shape of inc.double_conv.0.weight:", checkpoint['model_state_dict']['inc.double_conv.0.weight'].shape)
|
|
||||||
@@ -1,4 +0,0 @@
|
|||||||
import torch
|
|
||||||
from cloud_removal import DeepInpaintingStrategy
|
|
||||||
cloud_remover = DeepInpaintingStrategy()
|
|
||||||
print("Model channels:", list(cloud_remover.model.parameters())[0].shape[1])
|
|
||||||
@@ -1,10 +0,0 @@
|
|||||||
from cloud_removal import DeepInpaintingStrategy
|
|
||||||
import torch
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
cr = DeepInpaintingStrategy(model_path="cloud_removal_model/cloud_removal_unet_best.pth")
|
|
||||||
if cr.model is not None:
|
|
||||||
expected = list(cr.model.parameters())[0].shape[1]
|
|
||||||
print("Expected channels:", expected)
|
|
||||||
else:
|
|
||||||
print("Failed to load model")
|
|
||||||
@@ -1,53 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
import planetary_computer
|
|
||||||
import pystac_client
|
|
||||||
import odc.stac
|
|
||||||
import numpy as np
|
|
||||||
from shapely.geometry import Point, shape
|
|
||||||
from pyproj import Transformer
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.3, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-04-30"
|
|
||||||
catalog = pystac_client.Client.open(
|
|
||||||
"https://planetarycomputer.microsoft.com/api/stac/v1",
|
|
||||||
modifier=planetary_computer.sign_inplace,
|
|
||||||
)
|
|
||||||
items = list(catalog.search(collections=["sentinel-2-l2a"], bbox=bbox, datetime=time_range, query={"eo:cloud_cover": {"lt": 30}}).items())
|
|
||||||
items = sorted(items, key=lambda x: x.properties["eo:cloud_cover"])
|
|
||||||
|
|
||||||
gdf = gpd.read_file("train/ST_training_data_updated_1130points_new.shp")
|
|
||||||
gdf = gdf.to_crs("EPSG:32648")
|
|
||||||
|
|
||||||
# Find a point that fails. Let's just test a few points.
|
|
||||||
for idx, row in gdf.head(20).iterrows():
|
|
||||||
x_coord = row['geometry'].x
|
|
||||||
y_coord = row['geometry'].y
|
|
||||||
|
|
||||||
transformer = Transformer.from_crs("epsg:32648", "epsg:4326", always_xy=True)
|
|
||||||
lon, lat = transformer.transform(x_coord, y_coord)
|
|
||||||
point = Point(lon, lat)
|
|
||||||
|
|
||||||
filtered = []
|
|
||||||
for item in items:
|
|
||||||
if shape(item.geometry).contains(point):
|
|
||||||
filtered.append(item)
|
|
||||||
|
|
||||||
filtered = [planetary_computer.sign(item) for item in filtered]
|
|
||||||
if not filtered:
|
|
||||||
print(f"Point {idx}: NO ITEMS CONTAINS POINT!")
|
|
||||||
continue
|
|
||||||
|
|
||||||
ds = odc.stac.load(
|
|
||||||
filtered,
|
|
||||||
bands=["B02"],
|
|
||||||
x=(x_coord - 80, x_coord + 80),
|
|
||||||
y=(y_coord - 80, y_coord + 80),
|
|
||||||
crs="EPSG:32648",
|
|
||||||
resolution=10,
|
|
||||||
patch_url=planetary_computer.sign,
|
|
||||||
fail_on_error=False
|
|
||||||
).compute()
|
|
||||||
|
|
||||||
sums = ds["B02"].sum(dim=["x", "y"]).values
|
|
||||||
non_zero = (sums > 0).sum()
|
|
||||||
print(f"Point {idx}: {len(filtered)} items, {non_zero} non-zero time steps")
|
|
||||||
@@ -1,37 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
import planetary_computer
|
|
||||||
import pystac_client
|
|
||||||
import odc.stac
|
|
||||||
import sys
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.3, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-04-30"
|
|
||||||
catalog = pystac_client.Client.open("https://planetarycomputer.microsoft.com/api/stac/v1", modifier=planetary_computer.sign_inplace)
|
|
||||||
search = catalog.search(collections=["sentinel-2-l2a"], bbox=bbox, datetime=time_range, query={"eo:cloud_cover": {"lt": 30}})
|
|
||||||
items = list(search.items())
|
|
||||||
items = sorted(items, key=lambda x: x.properties.get("eo:cloud_cover", 100))[:4]
|
|
||||||
items = [planetary_computer.sign(item) for item in items]
|
|
||||||
|
|
||||||
gdf = gpd.read_file("train/ST_training_data_updated_1130points_new.shp")
|
|
||||||
gdf = gdf.to_crs("EPSG:32648")
|
|
||||||
|
|
||||||
row = gdf.iloc[0]
|
|
||||||
x, y_coord = row.geometry.x, row.geometry.y
|
|
||||||
point_bbox = [x - 80, y_coord - 80, x + 80, y_coord + 80]
|
|
||||||
|
|
||||||
patch_s2 = odc.stac.load(
|
|
||||||
items,
|
|
||||||
bands=["B02", "B03", "B04", "B08", "SCL"],
|
|
||||||
x=(x - 80, x + 80),
|
|
||||||
y=(y_coord - 80, y_coord + 80),
|
|
||||||
crs="EPSG:32648",
|
|
||||||
resolution=10,
|
|
||||||
patch_url=planetary_computer.sign,
|
|
||||||
fail_on_error=False
|
|
||||||
).compute()
|
|
||||||
|
|
||||||
print("patch_s2 vars:", patch_s2.data_vars)
|
|
||||||
if patch_s2.dims['x'] < 16 or patch_s2.dims['y'] < 16:
|
|
||||||
print("Too small:", patch_s2.dims)
|
|
||||||
else:
|
|
||||||
print("Success dimension:", patch_s2.dims)
|
|
||||||
@@ -1,3 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
gdf = gpd.read_file("train/ST_training_data_updated_1130points_new.shp")
|
|
||||||
print(gdf.head(1)['HT_code'])
|
|
||||||
@@ -1,3 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
gdf = gpd.read_file("train/ST_training_data_updated_1130points_new.shp")
|
|
||||||
print(gdf.columns)
|
|
||||||
@@ -1,23 +0,0 @@
|
|||||||
import numpy as np
|
|
||||||
from xgboost import XGBClassifier
|
|
||||||
from sklearn.datasets import make_classification
|
|
||||||
from sklearn.metrics import accuracy_score
|
|
||||||
|
|
||||||
print("🚀 Testing XGBoost with CUDA GPU...")
|
|
||||||
try:
|
|
||||||
X, y = make_classification(n_samples=10000, n_features=20, n_classes=2, random_state=42)
|
|
||||||
model = XGBClassifier(
|
|
||||||
n_estimators=100,
|
|
||||||
max_depth=10,
|
|
||||||
tree_method="hist",
|
|
||||||
device="cuda",
|
|
||||||
random_state=42,
|
|
||||||
verbosity=1
|
|
||||||
)
|
|
||||||
print("Training model...")
|
|
||||||
model.fit(X, y)
|
|
||||||
y_pred = model.predict(X)
|
|
||||||
acc = accuracy_score(y, y_pred)
|
|
||||||
print(f"✅ Training successful! Accuracy: {acc*100:.2f}%")
|
|
||||||
except Exception as e:
|
|
||||||
print(f"❌ Error during training: {e}")
|
|
||||||
@@ -1,3 +0,0 @@
|
|||||||
import torch
|
|
||||||
checkpoint = torch.load('cloud_removal_model/cloud_removal_unet_best.pth', map_location='cpu')
|
|
||||||
print(list(checkpoint['model_state_dict'].keys())[:5])
|
|
||||||
@@ -1,13 +0,0 @@
|
|||||||
import torch
|
|
||||||
from pathlib import Path
|
|
||||||
model = torch.load('cloud_removal_model/cloud_removal_unet_best.pth', map_location='cpu')
|
|
||||||
print(type(model))
|
|
||||||
print("hasattr inc:", hasattr(model, 'inc'))
|
|
||||||
if hasattr(model, 'inc'):
|
|
||||||
print("hasattr double_conv:", hasattr(model.inc, 'double_conv'))
|
|
||||||
if hasattr(model.inc, 'double_conv'):
|
|
||||||
print("in_channels:", model.inc.double_conv[0].in_channels)
|
|
||||||
else:
|
|
||||||
for name, param in model.named_parameters():
|
|
||||||
print(name, param.shape)
|
|
||||||
break
|
|
||||||
@@ -1,19 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
import planetary_computer
|
|
||||||
import pystac_client
|
|
||||||
import odc.stac
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.3, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-01-30"
|
|
||||||
catalog = pystac_client.Client.open("https://planetarycomputer.microsoft.com/api/stac/v1", modifier=planetary_computer.sign_inplace)
|
|
||||||
items = list(catalog.search(collections=["sentinel-2-l2a"], bbox=bbox, datetime=time_range).items())[:1]
|
|
||||||
items = [planetary_computer.sign(item) for item in items]
|
|
||||||
|
|
||||||
x = 561609
|
|
||||||
y = 1024183
|
|
||||||
try:
|
|
||||||
ds = odc.stac.load(items, bands=["B02"], crs="EPSG:32648", resolution=10, x=(x-80, x+80), y=(y-80, y+80))
|
|
||||||
print("Success with x/y:", ds.dims)
|
|
||||||
except Exception as e:
|
|
||||||
print("Error with x/y:", e)
|
|
||||||
|
|
||||||
@@ -1,18 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
import planetary_computer
|
|
||||||
import pystac_client
|
|
||||||
import odc.stac
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.3, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-04-30"
|
|
||||||
catalog = pystac_client.Client.open("https://planetarycomputer.microsoft.com/api/stac/v1", modifier=planetary_computer.sign_inplace)
|
|
||||||
# Get ALL items
|
|
||||||
items = list(catalog.search(collections=["sentinel-2-l2a"], bbox=bbox, datetime=time_range, query={"eo:cloud_cover": {"lt": 30}}).items())
|
|
||||||
items = [planetary_computer.sign(item) for item in items]
|
|
||||||
print(f"Total items: {len(items)}")
|
|
||||||
|
|
||||||
x = 561609
|
|
||||||
y = 1024183
|
|
||||||
ds = odc.stac.load(items, bands=["B02"], x=(x-80, x+80), y=(y-80, y+80), crs="EPSG:32648", resolution=10, patch_url=planetary_computer.sign).compute()
|
|
||||||
print("Time dimension size:", ds.dims['time'])
|
|
||||||
@@ -1,20 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
import planetary_computer
|
|
||||||
import pystac_client
|
|
||||||
import odc.stac
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.3, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-04-30"
|
|
||||||
catalog = pystac_client.Client.open("https://planetarycomputer.microsoft.com/api/stac/v1", modifier=planetary_computer.sign_inplace)
|
|
||||||
items = list(catalog.search(collections=["sentinel-2-l2a"], bbox=bbox, datetime=time_range, query={"eo:cloud_cover": {"lt": 30}}).items())
|
|
||||||
items = [planetary_computer.sign(item) for item in items]
|
|
||||||
|
|
||||||
x = 561609
|
|
||||||
y = 1024183
|
|
||||||
ds = odc.stac.load(items, bands=["B02"], x=(x-80, x+80), y=(y-80, y+80), crs="EPSG:32648", resolution=10, patch_url=planetary_computer.sign, fail_on_error=False).compute()
|
|
||||||
print("Original time size:", len(ds.time))
|
|
||||||
ds2 = ds.dropna(dim="time", how="all")
|
|
||||||
print("After dropna time size:", len(ds2.time))
|
|
||||||
print("B02 mean:", np.nanmean(ds2["B02"].values))
|
|
||||||
print("B02 non-nan count:", np.sum(~np.isnan(ds2["B02"].values)))
|
|
||||||
@@ -1,25 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
import planetary_computer
|
|
||||||
import pystac_client
|
|
||||||
import odc.stac
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.3, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-04-30"
|
|
||||||
catalog = pystac_client.Client.open("https://planetarycomputer.microsoft.com/api/stac/v1", modifier=planetary_computer.sign_inplace)
|
|
||||||
items = list(catalog.search(collections=["sentinel-2-l2a"], bbox=bbox, datetime=time_range, query={"eo:cloud_cover": {"lt": 30}}).items())
|
|
||||||
items = [planetary_computer.sign(item) for item in items]
|
|
||||||
|
|
||||||
x = 561609
|
|
||||||
y = 1024183
|
|
||||||
ds = odc.stac.load(items, bands=["B02", "B03", "B04", "B08", "SCL"], x=(x-80, x+80), y=(y-80, y+80), crs="EPSG:32648", resolution=10, patch_url=planetary_computer.sign, fail_on_error=False).compute()
|
|
||||||
|
|
||||||
print("Original shape:", ds["B02"].shape)
|
|
||||||
ds2 = ds.dropna(dim="time", how="all")
|
|
||||||
print("After dropna time size:", len(ds2.time))
|
|
||||||
if len(ds2.time) > 0:
|
|
||||||
ds2 = ds2.isel(time=slice(0, 4))
|
|
||||||
median = ds2["B02"].median(dim="time", skipna=True).values
|
|
||||||
print("Median shape:", median.shape)
|
|
||||||
print("Zeros in median:", np.sum(median == 0) / median.size)
|
|
||||||
print("NaNs in median:", np.sum(np.isnan(median)) / median.size)
|
|
||||||
@@ -1,15 +0,0 @@
|
|||||||
import sys
|
|
||||||
# Thêm đường dẫn hiện tại vào PYTHONPATH để import được new_import_ODC nếu cần
|
|
||||||
sys.path.append('.')
|
|
||||||
import warnings
|
|
||||||
warnings.filterwarnings('ignore')
|
|
||||||
from new_import_ODC import load_sen1
|
|
||||||
|
|
||||||
print("Testing load_sen1 with a short time range to speed up Dask compute...")
|
|
||||||
bbox = [105.5, 9.2, 106.4, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-01-31" # Short time range for fast testing
|
|
||||||
vh, vv = load_sen1(bbox, time_range)
|
|
||||||
print("VH shape:", vh.shape)
|
|
||||||
print("VV shape:", vv.shape)
|
|
||||||
print("VH CRS:", vh.rio.crs)
|
|
||||||
print("Success!")
|
|
||||||
@@ -1,49 +0,0 @@
|
|||||||
import warnings
|
|
||||||
warnings.filterwarnings('ignore')
|
|
||||||
|
|
||||||
def load_sen1(bbox, time_range):
|
|
||||||
import pystac_client
|
|
||||||
import planetary_computer
|
|
||||||
import odc.stac
|
|
||||||
|
|
||||||
catalog = pystac_client.Client.open(
|
|
||||||
"https://planetarycomputer.microsoft.com/api/stac/v1",
|
|
||||||
modifier=planetary_computer.sign_inplace,
|
|
||||||
)
|
|
||||||
|
|
||||||
search = catalog.search(
|
|
||||||
collections=["sentinel-1-rtc"],
|
|
||||||
bbox=bbox,
|
|
||||||
datetime=time_range,
|
|
||||||
)
|
|
||||||
items = list(search.items())
|
|
||||||
print("Found items:", len(items))
|
|
||||||
|
|
||||||
ds_s1 = odc.stac.load(
|
|
||||||
items,
|
|
||||||
bands=["vv", "vh"],
|
|
||||||
bbox=bbox,
|
|
||||||
crs="EPSG:32648",
|
|
||||||
resolution=10,
|
|
||||||
chunks={"x": 2048, "y": 2048, "time": 1}
|
|
||||||
)
|
|
||||||
|
|
||||||
ds_median = ds_s1.median(dim="time").compute()
|
|
||||||
vv = ds_median["vv"]
|
|
||||||
vh = ds_median["vh"]
|
|
||||||
|
|
||||||
vv = vv.expand_dims(dim="band")
|
|
||||||
vh = vh.expand_dims(dim="band")
|
|
||||||
|
|
||||||
vv = vv.rio.write_crs("EPSG:32648")
|
|
||||||
vh = vh.rio.write_crs("EPSG:32648")
|
|
||||||
|
|
||||||
return vh, vv
|
|
||||||
|
|
||||||
print("Testing load_sen1...")
|
|
||||||
bbox = [105.5, 9.2, 106.4, 10.0]
|
|
||||||
time_range = "2022-09-01/2023-10-01"
|
|
||||||
vh, vv = load_sen1(bbox, time_range)
|
|
||||||
print("VH shape:", vh.shape)
|
|
||||||
print("VV shape:", vv.shape)
|
|
||||||
print("Success!")
|
|
||||||
@@ -1,17 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
import planetary_computer
|
|
||||||
import pystac_client
|
|
||||||
import odc.stac
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.3, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-04-30"
|
|
||||||
catalog = pystac_client.Client.open("https://planetarycomputer.microsoft.com/api/stac/v1", modifier=planetary_computer.sign_inplace)
|
|
||||||
items = list(catalog.search(collections=["sentinel-2-l2a"], bbox=bbox, datetime=time_range, query={"eo:cloud_cover": {"lt": 30}}).items())[:4]
|
|
||||||
items = [planetary_computer.sign(item) for item in items]
|
|
||||||
|
|
||||||
x = 561609
|
|
||||||
y = 1024183
|
|
||||||
ds = odc.stac.load(items, bands=["B02", "B03", "B04", "B08"], x=(x-80, x+80), y=(y-80, y+80), crs="EPSG:32648", resolution=10, patch_url=planetary_computer.sign).compute()
|
|
||||||
print("B04 nanmean:", np.nanmean(ds["B04"].values))
|
|
||||||
print("B04 nanmax:", np.nanmax(ds["B04"].values))
|
|
||||||
@@ -1,41 +0,0 @@
|
|||||||
import new_import_ODC
|
|
||||||
importlib = __import__('importlib')
|
|
||||||
importlib.reload(new_import_ODC)
|
|
||||||
from new_import_ODC import *
|
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
date_range = ("2022-09-01", "2022-10-01")
|
|
||||||
longtitude_range = (105.86, 105.94)
|
|
||||||
latitude_range = (9.65, 9.69)
|
|
||||||
coordinates = (longtitude_range, latitude_range)
|
|
||||||
|
|
||||||
print("Loading S2...")
|
|
||||||
data = load_data(None, date_range, longtitude_range, latitude_range)
|
|
||||||
result = mask_clean(data)
|
|
||||||
ds1 = calculate_indices(result, index="NDVI", satellite_mission="s2")
|
|
||||||
ndvi = ds1["NDVI"]
|
|
||||||
time_split = [
|
|
||||||
slice("2022-09-01", "2023-01-01"),
|
|
||||||
slice("2023-01-01", "2023-05-01"),
|
|
||||||
slice("2023-05-01", "2023-07-01"),
|
|
||||||
slice("2023-07-01", "2022-10-01"),
|
|
||||||
]
|
|
||||||
fill_nan_ndvi = fill_nan(ndvi, time_split)
|
|
||||||
average_ndvi = fill_nan_ndvi.resample(time="1M").mean().compute()
|
|
||||||
|
|
||||||
print("Loading S1...")
|
|
||||||
dsvh, dsvv = load_data_sen1(None, date_range, coordinates)
|
|
||||||
average_vv = calculate_average(dsvv, time_pattern='1M')
|
|
||||||
average_vh = calculate_average(dsvh, time_pattern='1M')
|
|
||||||
|
|
||||||
train = load_train_data("train/ST_training_data_updated_1130points_new.shp")
|
|
||||||
point = train.iloc[0]
|
|
||||||
|
|
||||||
ndvi_val = average_ndvi.sel(x=point.geometry.x, y=point.geometry.y, method='nearest').values
|
|
||||||
vh_val = average_vh.sel(x=point.geometry.x, y=point.geometry.y, method='nearest').values
|
|
||||||
vv_val = average_vv.sel(x=point.geometry.x, y=point.geometry.y, method='nearest').values
|
|
||||||
|
|
||||||
print("NDVI shape:", ndvi_val.shape, "ndim:", ndvi_val.ndim)
|
|
||||||
print("VH shape:", vh_val.shape, "ndim:", vh_val.ndim)
|
|
||||||
print("VV shape:", vv_val.shape, "ndim:", vv_val.ndim)
|
|
||||||
|
|
||||||
@@ -1,42 +0,0 @@
|
|||||||
import geopandas as gpd
|
|
||||||
import planetary_computer
|
|
||||||
import pystac_client
|
|
||||||
import odc.stac
|
|
||||||
import numpy as np
|
|
||||||
import time
|
|
||||||
from shapely.geometry import Point, box, shape
|
|
||||||
|
|
||||||
bbox = [105.5, 9.2, 106.3, 10.0]
|
|
||||||
time_range = "2023-01-01/2023-04-30"
|
|
||||||
catalog = pystac_client.Client.open("https://planetarycomputer.microsoft.com/api/stac/v1", modifier=planetary_computer.sign_inplace)
|
|
||||||
items = list(catalog.search(collections=["sentinel-2-l2a"], bbox=bbox, datetime=time_range, query={"eo:cloud_cover": {"lt": 30}}).items())
|
|
||||||
|
|
||||||
x = 561609
|
|
||||||
y = 1024183
|
|
||||||
|
|
||||||
start = time.time()
|
|
||||||
# Filter items by spatial intersection
|
|
||||||
from pyproj import Transformer
|
|
||||||
# The items geometry are in EPSG:4326 (lon, lat)
|
|
||||||
# Our x, y are in EPSG:32648
|
|
||||||
transformer = Transformer.from_crs("epsg:32648", "epsg:4326", always_xy=True)
|
|
||||||
lon, lat = transformer.transform(x, y)
|
|
||||||
point = Point(lon, lat)
|
|
||||||
|
|
||||||
filtered_items = []
|
|
||||||
for item in items:
|
|
||||||
geom = shape(item.geometry)
|
|
||||||
if geom.contains(point):
|
|
||||||
filtered_items.append(item)
|
|
||||||
|
|
||||||
filtered_items = sorted(filtered_items, key=lambda x: x.properties["eo:cloud_cover"])
|
|
||||||
|
|
||||||
print("Original items:", len(items))
|
|
||||||
print("Filtered items:", len(filtered_items))
|
|
||||||
print("Time to filter:", time.time() - start)
|
|
||||||
|
|
||||||
start = time.time()
|
|
||||||
filtered_items = [planetary_computer.sign(item) for item in filtered_items]
|
|
||||||
ds = odc.stac.load(filtered_items[:4], bands=["B02"], x=(x-80, x+80), y=(y-80, y+80), crs="EPSG:32648", resolution=10, patch_url=planetary_computer.sign, fail_on_error=False).compute()
|
|
||||||
print("Time to load 4 items:", time.time() - start)
|
|
||||||
print(ds["B02"].shape)
|
|
||||||
@@ -1,5 +0,0 @@
|
|||||||
from train_cloud_removal import UNet
|
|
||||||
model = UNet(in_channels=6, out_channels=4)
|
|
||||||
print(hasattr(model, 'inc'))
|
|
||||||
print(hasattr(model, 'conv1'))
|
|
||||||
print(list(model.parameters())[0].shape)
|
|
||||||
@@ -1,13 +0,0 @@
|
|||||||
import numpy as np
|
|
||||||
|
|
||||||
y = []
|
|
||||||
# simulate appending 1130 labels
|
|
||||||
for i in range(1130):
|
|
||||||
y.append(i % 5)
|
|
||||||
|
|
||||||
y = np.array(y)
|
|
||||||
unique_labels = sorted(list(np.unique(y)))
|
|
||||||
label_map = {lbl: i for i, lbl in enumerate(unique_labels)}
|
|
||||||
y_mapped = np.array([label_map[l] for l in y])
|
|
||||||
|
|
||||||
print(len(y), len(y_mapped))
|
|
||||||
+1135
File diff suppressed because it is too large
Load Diff
File diff suppressed because one or more lines are too long
+102
-102
File diff suppressed because one or more lines are too long
File diff suppressed because one or more lines are too long
@@ -0,0 +1,638 @@
|
|||||||
|
# GEMINI PROJECT CONTEXT - Land Classification & Remote Sensing System
|
||||||
|
|
||||||
|
**Last Updated**: March 26, 2026
|
||||||
|
**Project Location**: `/home/x79/remote-sensing`
|
||||||
|
**Purpose**: Complete land classification and environmental monitoring system using satellite remote sensing for Vietnam
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 📋 PROJECT OVERVIEW
|
||||||
|
|
||||||
|
### High-Level Purpose & Problem Domain
|
||||||
|
- **Core Task**: Classify land use/land cover (8 land classes) in Vietnam using multispectral Sentinel-2 and radar Sentinel-1 data from Microsoft Planetary Computer
|
||||||
|
- **Geographic Focus**: Vietnam provinces/regions with bounding-box (bbox) based Area-of-Interest (AOI) selection
|
||||||
|
- **Key Capabilities**:
|
||||||
|
- Dynamic training with user-selected regions and time periods
|
||||||
|
- Pixel-wise inference (prediction) on new regions
|
||||||
|
- Cloud removal using 7 different strategies
|
||||||
|
- NDVI time-series forecasting and change detection workflows
|
||||||
|
- Auto-generated HTML reports with visualizations
|
||||||
|
- Batch processing of multiple regions
|
||||||
|
- Model lifecycle management (save, load, validate, delete)
|
||||||
|
|
||||||
|
### Data Pipeline
|
||||||
|
```
|
||||||
|
Sentinel-2 (optical) + Sentinel-1 (SAR)
|
||||||
|
↓
|
||||||
|
[Feature Extraction: 4 modes - simple (3) / temporal (39) / extended (15) / odc (8)]
|
||||||
|
↓
|
||||||
|
[Model Training: XGBoost, RF, SVM, CNN, Swin-UNet, MobileNet-LRASPP]
|
||||||
|
↓
|
||||||
|
[Prediction: Pixel-wise classification]
|
||||||
|
↓
|
||||||
|
[Output: GeoTIFF + PNG preview + HTML report + JSON metadata]
|
||||||
|
```
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 🏗️ SYSTEM ARCHITECTURE
|
||||||
|
|
||||||
|
### Core Technology Stack
|
||||||
|
- **Backend**: FastAPI (~4200 lines in `api_server.py`)
|
||||||
|
- **ML Training**: scikit-learn (XGBoost, RF, SVM, DT) + PyTorch (CNN, Swin-UNet, MobileNet)
|
||||||
|
- **Geospatial**: rasterio, rioxarray, geopandas, xarray, odc.stac
|
||||||
|
- **Data Access**: Microsoft Planetary Computer STAC API (Sentinel-2 L2A, Sentinel-1 RTC)
|
||||||
|
- **Frontend**: HTML + Leaflet.js (map drawing) + Fetch API + Chart.js
|
||||||
|
- **GPU Support**: PyTorch with CUDA 12.x (optional fallback to CPU)
|
||||||
|
|
||||||
|
### Folder Structure
|
||||||
|
```
|
||||||
|
remote-sensing/
|
||||||
|
├── Core Backend
|
||||||
|
│ ├── api_server.py # FastAPI app (~4200 LOC, 70+ endpoints)
|
||||||
|
│ ├── train_module.py # Training pipeline engine
|
||||||
|
│ ├── feature_extractor.py # Unified feature extraction (4 modes)
|
||||||
|
│ ├── model_manager.py # Model lifecycle management
|
||||||
|
│ ├── cloud_removal.py # 7 cloud removal strategies
|
||||||
|
│ ├── report_generator.py # Auto HTML/PNG report generation
|
||||||
|
│ ├── generate_previews.py # GeoTIFF → PNG conversion
|
||||||
|
│ │
|
||||||
|
├── Utilities & Lookup
|
||||||
|
│ ├── vietnam_provinces.py # Province bboxes & metadata
|
||||||
|
│ ├── vietnam_provinces_merged.py # 32-province variant
|
||||||
|
│ ├── utils.py # Geospatial helper functions
|
||||||
|
│ ├── create_odc_metadata.py # Metadata generator utility
|
||||||
|
│ │
|
||||||
|
├── Frontend Pages (HTML)
|
||||||
|
│ ├── index.html # Main dashboard hub
|
||||||
|
│ ├── training_interface.html # Training UI
|
||||||
|
│ ├── prediction_interface.html # Prediction UI
|
||||||
|
│ ├── batch_interface.html # Batch processing UI
|
||||||
|
│ ├── ndvi_interface.html # NDVI time-series UI
|
||||||
|
│ ├── dashboard.html # Analytics dashboard
|
||||||
|
│ ├── reports_interface.html # Reports management
|
||||||
|
│ ├── change_detection_interface.html # Change detection UI
|
||||||
|
│ ├── cloud_training_interface.html # Cloud removal training UI
|
||||||
|
│ │
|
||||||
|
├── Tests & Notebooks
|
||||||
|
│ ├── test_*.py # Unit & integration tests
|
||||||
|
│ ├── 01.train_ODC*.ipynb # Training notebooks
|
||||||
|
│ ├── 02.predict_ODC.ipynb # Prediction notebooks
|
||||||
|
│ ├── cloud_removal_train.ipynb # Cloud removal training
|
||||||
|
│ │
|
||||||
|
├── Model Storage & Caches
|
||||||
|
│ ├── model_train/ # Trained models (*.joblib, *.pth)
|
||||||
|
│ │ ├── model_odc.joblib # Legacy GridSearchCV model
|
||||||
|
│ │ ├── model_*_info.json # Metadata sidecar files
|
||||||
|
│ ├── cloud_removal_model/ # Cloud removal U-Net models (.pth)
|
||||||
|
│ ├── predictions/ # Prediction output (GeoTIFF + PNG)
|
||||||
|
│ ├── reports/ # Generated HTML reports
|
||||||
|
│ ├── dataset_cache/ # Cached Sentinel data (optional)
|
||||||
|
│ │
|
||||||
|
├── Config & Documentation
|
||||||
|
│ ├── requirement.txt # Python dependencies
|
||||||
|
│ ├── requirements_api.txt # API-specific deps
|
||||||
|
│ ├── IMPLEMENTATION_SUMMARY.md # Model manager summary
|
||||||
|
│ ├── MODEL_MANAGER_GUIDE.md # Full model management guide
|
||||||
|
│ ├── NDVI_FORECAST_METHODOLOGY.md # NDVI algorithm docs
|
||||||
|
│ ├── CLOUD_TRAINING_GUIDE.md # Cloud removal training guide
|
||||||
|
│ └── [Other guides & docs]
|
||||||
|
```
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 🔧 MAIN MODULES & RESPONSIBILITIES
|
||||||
|
|
||||||
|
| **Module** | **File(s)** | **Key Responsibility** |
|
||||||
|
|---|---|---|
|
||||||
|
| **API Server** | `api_server.py` | FastAPI app with 70+ endpoints; routes all training, prediction, batch, cloud removal, dashboard, reports, model management tasks |
|
||||||
|
| **Training Engine** | `train_module.py` | Complete training pipeline: fetch data → feature extraction → train/test split → model training → evaluation → save with metadata |
|
||||||
|
| **Feature Extraction** | `feature_extractor.py` | Standardized feature extraction with 4 modes: simple, temporal, extended, odc; used by both training and prediction |
|
||||||
|
| **Model Manager** | `model_manager.py` | Lifecycle management: list, load, save, validate, delete models; handles metadata JSON; auto-detects CNN/PyTorch models |
|
||||||
|
| **Cloud Removal** | `cloud_removal.py` | 7 cloud removal strategies: classic (3-step), temporal_only, median_composite, none, speckle filter, ML inpainting, deep learning U-Net |
|
||||||
|
| **Report Generator** | `report_generator.py` | Auto-generates HTML/PNG reports with confusion matrices, class distributions, accuracy trends |
|
||||||
|
| **Preview Generator** | `generate_previews.py` | Converts GeoTIFF outputs to PNG previews (NDVI or classification rasters) |
|
||||||
|
| **Province Lookup** | `vietnam_provinces*.py` | Lookup tables for 32+ Vietnamese provinces with bboxes and region grouping |
|
||||||
|
| **Utilities** | `utils.py` | Geospatial helper functions (load GeoDataFrames, etc.) |
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 📊 END-TO-END WORKFLOWS
|
||||||
|
|
||||||
|
### 1. TRAINING WORKFLOW
|
||||||
|
```
|
||||||
|
User Input → Training Configuration
|
||||||
|
↓
|
||||||
|
API Endpoint: POST /api/training/start
|
||||||
|
↓
|
||||||
|
train_module.py: train_model()
|
||||||
|
1. Fetch Sentinel-2 & Sentinel-1 from Planetary Computer STAC
|
||||||
|
2. Apply cloud mask (SCL band: clouds, shadows, cirrus masked)
|
||||||
|
3. Extract features via FeatureExtractor (mode: simple/temporal/extended/odc)
|
||||||
|
4. Train/test split (default 0.2)
|
||||||
|
5. Train selected model type (XGBoost, RF, CNN, Swin-UNet, MobileNet)
|
||||||
|
6. Evaluate: accuracy, precision, recall, F1, confusion matrix
|
||||||
|
↓
|
||||||
|
model_manager.py: Save model + JSON metadata
|
||||||
|
↓
|
||||||
|
report_generator.py: Auto-generate HTML training report
|
||||||
|
↓
|
||||||
|
Return: {model_filename, accuracy_metrics, training_time}
|
||||||
|
```
|
||||||
|
|
||||||
|
**Key Metadata Saved**:
|
||||||
|
```json
|
||||||
|
{
|
||||||
|
"timestamp": "2026-03-26T14:30:00",
|
||||||
|
"model_type": "xgboost",
|
||||||
|
"feature_mode": "temporal",
|
||||||
|
"n_features": 39,
|
||||||
|
"n_classes": 8,
|
||||||
|
"features": ["NDVI_t1", "NDVI_t2", ..., "NDWI_t1", ...],
|
||||||
|
"test_accuracy": 0.85,
|
||||||
|
"train_accuracy": 0.92,
|
||||||
|
"bbox": [105.6, 9.3, 106.2, 9.8],
|
||||||
|
"time_range": "2023-03-01/2023-05-31",
|
||||||
|
"resolution": 20,
|
||||||
|
"data_source": "Microsoft Planetary Computer STAC"
|
||||||
|
}
|
||||||
|
```
|
||||||
|
|
||||||
|
### 2. PREDICTION WORKFLOW
|
||||||
|
```
|
||||||
|
User Input → Prediction Configuration (model_filename, bbox, time_range, cloud_strategy)
|
||||||
|
↓
|
||||||
|
API Endpoint: POST /api/predict or POST /api/predict/with-ndvi
|
||||||
|
↓
|
||||||
|
run_prediction() function:
|
||||||
|
1. Load model via model_manager.py (retrieves metadata, feature requirements)
|
||||||
|
2. Fetch Sentinel-2 & Sentinel-1 for new region
|
||||||
|
3. Apply chosen cloud_removal_method (classic/temporal_only/median_composite/none/deep_learning)
|
||||||
|
4. Extract features matching model's metadata requirements
|
||||||
|
5. Auto-adjust if feature count mismatch (pad/trim)
|
||||||
|
6. Predict land class for each pixel
|
||||||
|
7. (Optional) Calculate NDVI: (NIR - Red) / (NIR + Red)
|
||||||
|
8. Save outputs: GeoTIFF + PNG preview
|
||||||
|
↓
|
||||||
|
generate_previews.py: Create PNG from GeoTIFF
|
||||||
|
↓
|
||||||
|
report_generator.py: Generate prediction report
|
||||||
|
↓
|
||||||
|
Return: {prediction_file, ndvi_file, class_distribution, statistics}
|
||||||
|
```
|
||||||
|
|
||||||
|
### 3. BATCH PROCESSING WORKFLOW
|
||||||
|
```
|
||||||
|
User uploads CSV with multiple regions:
|
||||||
|
(name, min_lon, min_lat, max_lon, max_lat, start_date, end_date, max_scenes, cloud_cover, resolution)
|
||||||
|
↓
|
||||||
|
API Endpoint: POST /api/batch/start
|
||||||
|
↓
|
||||||
|
Enqueue all regions; process sequentially
|
||||||
|
↓
|
||||||
|
For each region: Run same prediction workflow
|
||||||
|
↓
|
||||||
|
Track status per region: Queued → Running → Completed/Failed
|
||||||
|
↓
|
||||||
|
UI shows progress bar, auto-retry on failure (max 3 retries)
|
||||||
|
↓
|
||||||
|
Return: Bulk results with per-region status & output files
|
||||||
|
```
|
||||||
|
|
||||||
|
### 4. CLOUD REMOVAL WORKFLOW
|
||||||
|
```
|
||||||
|
User selects cloud_removal_method in prediction config:
|
||||||
|
↓
|
||||||
|
cloud_removal.py: process_cloud_removal()
|
||||||
|
↓
|
||||||
|
Strategy Selection:
|
||||||
|
• 'classic': temporal interpolation → median composite → spatial interpolation (3-step)
|
||||||
|
• 'temporal_only': ffill + bfill across time dimension (fast, good for many scenes)
|
||||||
|
• 'median_composite': Prioritize median across scenes (best for noise reduction)
|
||||||
|
• 'none': Keep original, just fill NaN with 0
|
||||||
|
• 'deep': Use trained U-Net model (S2 cloudy + S1 → clean S2)
|
||||||
|
• 'ml_inpainting': KNN or Random Forest based inpainting
|
||||||
|
• 'speckle_filter': Reduce radar noise
|
||||||
|
↓
|
||||||
|
Return cleaned Sentinel-2 data for subsequent feature extraction
|
||||||
|
```
|
||||||
|
|
||||||
|
### 5. NDVI TIME-SERIES WORKFLOW
|
||||||
|
```
|
||||||
|
User requests NDVI calculation (bbox + time_range + aggregation)
|
||||||
|
↓
|
||||||
|
API Endpoint: POST /api/ndvi/timeseries or /api/ndvi/predict-timeseries
|
||||||
|
↓
|
||||||
|
Load Sentinel-2 (B04 Red, B08 NIR)
|
||||||
|
↓
|
||||||
|
Calculate NDVI = (NIR - Red) / (NIR + Red + 0.00001)
|
||||||
|
↓
|
||||||
|
Resample to monthly or user-defined aggregation
|
||||||
|
↓
|
||||||
|
Export as GeoTIFF + PNG visualization
|
||||||
|
↓
|
||||||
|
Show time-series graph & statistics (mean, min, max, std, trend)
|
||||||
|
```
|
||||||
|
|
||||||
|
### 6. CHANGE DETECTION WORKFLOW
|
||||||
|
```
|
||||||
|
User selects: model + current_period + prediction_period
|
||||||
|
↓
|
||||||
|
API Endpoint: POST /api/change-detection/compare-periods
|
||||||
|
↓
|
||||||
|
Run prediction for both time periods
|
||||||
|
↓
|
||||||
|
Compute difference map (current - prediction)
|
||||||
|
↓
|
||||||
|
Classify changes: increased vegetation, decreased vegetation, stable
|
||||||
|
↓
|
||||||
|
Generate change map GeoTIFF + report with statistics
|
||||||
|
```
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 🌐 API ENDPOINTS SUMMARY (70+ endpoints)
|
||||||
|
|
||||||
|
### Model Management
|
||||||
|
- `GET /api/models/list` - List all trained models with metadata
|
||||||
|
- `GET /api/models/{filename}/info` - Get model details
|
||||||
|
- `GET /api/models/{filename}/validate` - Validate model integrity
|
||||||
|
- `DELETE /api/models/{filename}` - Delete model file
|
||||||
|
|
||||||
|
### Training APIs
|
||||||
|
- `POST /api/training/start` - Start land classification training
|
||||||
|
- `GET /api/training/status` - Get training progress
|
||||||
|
- `POST /api/training/stop` - Cancel ongoing training
|
||||||
|
- `POST /api/cloud-removal/train` - Train cloud removal U-Net
|
||||||
|
|
||||||
|
### Prediction APIs
|
||||||
|
- `POST /api/predict` - Standard prediction (classification only)
|
||||||
|
- `POST /api/predict/with-ndvi` - Prediction with NDVI export
|
||||||
|
- `POST /api/change-detection/compare-periods` - Change detection
|
||||||
|
- `GET /api/prediction/status` - Check prediction progress
|
||||||
|
- `GET /api/predictions/list` - List prediction outputs
|
||||||
|
- `GET /api/predictions/download/{filename}` - Download prediction file
|
||||||
|
- `GET /api/predictions/preview/{filename}` - View PNG preview
|
||||||
|
|
||||||
|
### Batch Processing
|
||||||
|
- `POST /api/batch/start` - Enqueue multiple predictions from CSV
|
||||||
|
- `GET /api/batch/status` - Check batch queue
|
||||||
|
- `GET /api/batch/results/{batch_id}` - Retrieve batch results
|
||||||
|
- `POST /api/batch/cancel/{batch_id}` - Cancel batch job
|
||||||
|
|
||||||
|
### Cloud Removal
|
||||||
|
- `GET /api/cloud-removal/methods` - List available strategies
|
||||||
|
- `GET /api/cloud-removal/models` - List trained .pth models
|
||||||
|
- `POST /api/cloud-removal/upload` - Upload .pth cloud removal model
|
||||||
|
- `DELETE /api/cloud-removal/models/{filename}` - Delete cloud removal model
|
||||||
|
|
||||||
|
### Dashboard & Reports
|
||||||
|
- `GET /api/dashboard/statistics` - Overall system stats
|
||||||
|
- `GET /api/dashboard/accuracy-trends` - Accuracy over time
|
||||||
|
- `GET /api/dashboard/class-distribution/{model_filename}` - Class distribution
|
||||||
|
- `GET /api/reports/list` - List generated reports
|
||||||
|
- `GET /api/reports/view/{filename}` - View HTML report
|
||||||
|
- `GET /api/reports/download/{filename}` - Download report
|
||||||
|
- `DELETE /api/reports/delete/{filename}` - Delete report
|
||||||
|
|
||||||
|
### Provinces & Utilities
|
||||||
|
- `GET /api/provinces/list` - List all Vietnamese provinces
|
||||||
|
- `GET /api/provinces/by-region` - Group provinces by region
|
||||||
|
- `GET /api/provinces/{province_name}/bbox` - Get province bbox
|
||||||
|
- `GET /api/provinces/search/{query}` - Search province by name
|
||||||
|
- `GET /api/provinces-32/*` - Alternative 32-province variant
|
||||||
|
- `GET /api/network/check` - Check connectivity to Planetary Computer
|
||||||
|
- `GET /api/cache/info` - Show cache statistics
|
||||||
|
- `POST /api/cache/clear` - Clear local cache
|
||||||
|
|
||||||
|
### NDVI & Time-Series
|
||||||
|
- `POST /api/ndvi/timeseries` - Calculate NDVI time-series
|
||||||
|
- `POST /api/ndvi/predict-timeseries` - NDVI prediction/forecast
|
||||||
|
- `POST /api/ndvi/forecast` - NDVI forecasting
|
||||||
|
|
||||||
|
### File Management
|
||||||
|
- `GET /api/training/files` - List training files
|
||||||
|
- `GET /api/overlay/shapefiles` - List available shapefiles
|
||||||
|
- `GET /api/training/shapefile/{filename}/labels` - Get shapefile labels
|
||||||
|
- `POST /api/land-classification/upload` - Upload custom model
|
||||||
|
- `POST /api/cloud-removal/upload` - Upload cloud removal model
|
||||||
|
|
||||||
|
### Frontend Routes (Serve HTML)
|
||||||
|
- `GET /` - Main dashboard
|
||||||
|
- `GET /training` - Training interface
|
||||||
|
- `GET /prediction` - Prediction interface
|
||||||
|
- `GET /dashboard` - Analytics dashboard
|
||||||
|
- `GET /batch` - Batch processing UI
|
||||||
|
- `GET /ndvi` - NDVI time-series UI
|
||||||
|
- `GET /reports` - Reports management
|
||||||
|
- `GET /cloud-training` - Cloud removal training
|
||||||
|
- `GET /change-detection` - Change detection UI
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 💾 DATA INPUTS / OUTPUTS & FOLDER CONVENTIONS
|
||||||
|
|
||||||
|
### Input Data Sources
|
||||||
|
- **Sentinel-2 L2A** from Microsoft Planetary Computer STAC API
|
||||||
|
- Bands: B02 (blue), B03 (green), B04 (red), B08 (NIR), B11 (SWIR), SCL (cloud mask)
|
||||||
|
- Resolution: 10m or 20m (user selectable)
|
||||||
|
- Collection: `sentinel-2-l2a`
|
||||||
|
|
||||||
|
- **Sentinel-1 RTC** from Planetary Computer
|
||||||
|
- Bands: VH, VV (radar polarizations)
|
||||||
|
- Converted to dB scale: `10 * log10(intensity)`
|
||||||
|
- Collection: `sentinel-1-rtc`
|
||||||
|
|
||||||
|
- **Training Labels**: User-provided shapefiles with pixel-level class labels
|
||||||
|
|
||||||
|
### Output File Structure
|
||||||
|
```
|
||||||
|
predictions/
|
||||||
|
├── prediction_YYYYMMDD_HHMMSS.tif # Classification GeoTIFF
|
||||||
|
├── prediction_YYYYMMDD_HHMMSS.png # PNG preview
|
||||||
|
├── ndvi_YYYYMMDD_HHMMSS.tif # NDVI raster
|
||||||
|
├── ndvi_YYYYMMDD_HHMMSS.png # NDVI preview
|
||||||
|
|
||||||
|
reports/
|
||||||
|
├── training_report_*.html # Auto training reports
|
||||||
|
├── prediction_report_*.html # Auto prediction reports
|
||||||
|
|
||||||
|
model_train/
|
||||||
|
├── model_odc.joblib # Legacy model
|
||||||
|
├── model_odc_info.json # Metadata
|
||||||
|
├── model_xgboost_*.joblib # XGBoost models
|
||||||
|
├── model_xgboost_*_info.json # Metadata
|
||||||
|
├── model_cnn_*.joblib # CNN models
|
||||||
|
├── model_cnn_*_info.json # Metadata
|
||||||
|
|
||||||
|
cloud_removal_model/
|
||||||
|
├── cloud_removal_unet_best.pth # Trained U-Net
|
||||||
|
├── *.pth # Custom models
|
||||||
|
├── *.json # Model metadata
|
||||||
|
```
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 🔌 EXTERNAL DEPENDENCIES & PLATFORMS
|
||||||
|
|
||||||
|
### Critical External Services
|
||||||
|
- **Microsoft Planetary Computer** (STAC API)
|
||||||
|
- Hosts Sentinel-2 L2A and Sentinel-1 RTC archives
|
||||||
|
- URL: `https://planetarycomputer.microsoft.com/api/stac/v1`
|
||||||
|
- Auto-signed access tokens via `planetary_computer.sign_inplace`
|
||||||
|
- Network connectivity check: `GET /api/network/check`
|
||||||
|
|
||||||
|
### Key Python Libraries
|
||||||
|
- **Geospatial**: rasterio, rioxarray, geopandas, shapely, Cartopy, folium, ipyleaflet
|
||||||
|
- **Data Processing**: numpy, pandas, xarray, dask
|
||||||
|
- **ML**: scikit-learn, xgboost
|
||||||
|
- **Deep Learning**: torch, torchvision
|
||||||
|
- **Web**: fastapi, uvicorn, pydantic
|
||||||
|
- **Visualization**: matplotlib, Pillow (PIL)
|
||||||
|
- **Document Gen**: markdown, Pillow
|
||||||
|
|
||||||
|
### GPU Support
|
||||||
|
- PyTorch with CUDA 12.x (optional; falls back to CPU)
|
||||||
|
- Benefits Swin-UNet and CNN models (10-100x speedup)
|
||||||
|
- CPU training for XGBoost/RF typically <1 hour; deep models need GPU for reasonable speed
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## ⚙️ FEATURE EXTRACTION MODES (CRITICAL)
|
||||||
|
|
||||||
|
Train and prediction **MUST** use same feature mode and dimension; metadata auto-detects this.
|
||||||
|
|
||||||
|
| Mode | # Features | Description | Best For | Training Time |
|
||||||
|
|---|---|---|---|---|
|
||||||
|
| **simple** | 3 | NDVI_mean, VH_db_mean, VV_db_mean | Fast iteration, baseline | ~5-10 min |
|
||||||
|
| **temporal** | 39 | NDVI/NDWI/NDBI across 13 months + radar stats | High accuracy (~85%+) | ~30-60 min |
|
||||||
|
| **extended** | 15 | NDVI/NDWI/NDBI stats (mean/std/min/max) + radar | Balanced speed/accuracy | ~15-30 min |
|
||||||
|
| **odc** | 8 | NDVI stats + NDWI/NDBI/EVI mean (legacy ODC mode) | Legacy compatibility | ~10-20 min |
|
||||||
|
|
||||||
|
**Critical**: If feature mode = "temporal" (39 features) at training, prediction MUST extract 39 features. System auto-detects from metadata but will fail if mismatched.
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 🎯 OPERATIONAL NOTES & CONSTRAINTS
|
||||||
|
|
||||||
|
### Performance Limits
|
||||||
|
1. **Planetary Computer Timeout Issues**
|
||||||
|
- Large bbox (>10km × 10km) + long time range (>1 month) + high max_scenes → timeouts
|
||||||
|
- **Solution**: Progressive loading (subdivide bbox), reduce time window, reduce max_scenes
|
||||||
|
- **Safe Settings**: bbox ≤ 10km × 10km, time ≤ 1 month, max_scenes ≤ 12
|
||||||
|
|
||||||
|
2. **Memory Usage**
|
||||||
|
- Temporal mode (39 features) requires ~2-3x RAM vs simple mode
|
||||||
|
- Large regions: reduce resolution (10m → 20m) or split into sub-tiles
|
||||||
|
- Batch processing: sequential (one region at a time due to API limits)
|
||||||
|
|
||||||
|
3. **GPU Training**
|
||||||
|
- Swin-UNet: ~15-60 min on GPU vs ~2-4 hours on CPU
|
||||||
|
- CNN: ~10-30 min on GPU vs ~1-2 hours on CPU
|
||||||
|
- XGBoost/RF: CPU-bound; GPU not beneficial
|
||||||
|
|
||||||
|
### Data Quality Issues
|
||||||
|
1. **Cloud Cover**
|
||||||
|
- SCL band values: 3=cloud shadow, 8=cloud medium, 9=cloud high, 10=cirrus
|
||||||
|
- Recommend multiple scenes (≥5) for temporal aggregation
|
||||||
|
- Cloud removal strategy critical—test different approaches
|
||||||
|
|
||||||
|
2. **Radar Data (Sentinel-1)**
|
||||||
|
- Not always available for all regions/dates
|
||||||
|
- System gracefully falls back to zeros if unavailable
|
||||||
|
- Safe for "extended" & "odc" modes that have radar fallback
|
||||||
|
|
||||||
|
3. **Feature Mode Mismatch**
|
||||||
|
- Model trained with "temporal" (39 features) needs 39-dim input
|
||||||
|
- System auto-adjusts (pads/trims) from metadata but may degrade accuracy
|
||||||
|
- **Best Practice**: Align feature mode explicitly; don't mix
|
||||||
|
|
||||||
|
### Known Caveats
|
||||||
|
1. **Legacy Model (model_odc.joblib)**: Hardcoded 39 temporal features; auto-detected via `model_odc_info.json`
|
||||||
|
2. **Metadata Consistency**: Old models may lack `.json` sidecar; system generates default (may be incorrect)
|
||||||
|
3. **Batch Processing**: Sequential only; large batches (100+ regions) take hours
|
||||||
|
4. **Change Detection**: Simple differencing approach; requires same model & feature mode for both periods
|
||||||
|
5. **Rate Limiting**: Planetary Computer may rate-limit if too many concurrent requests
|
||||||
|
|
||||||
|
### Recommended Best Practices
|
||||||
|
- Test model on small bbox first (2km × 2km, 1 week, 3 scenes)
|
||||||
|
- Use "simple" mode for fast iteration, "temporal" for best accuracy (85%+)
|
||||||
|
- Store metadata JSON alongside model file (sidecar pattern)
|
||||||
|
- Version control: record feature_mode & n_features in every training
|
||||||
|
- Monitor training accuracy; retrain if <70% accuracy
|
||||||
|
- Cache Sentinel data locally to avoid repeated downloads
|
||||||
|
- Use "median_composite" cloud strategy if >5 scenes; "temporal_only" if 3-4 scenes
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 🔍 TEST COVERAGE MAP
|
||||||
|
|
||||||
|
| Test File | Coverage | Status |
|
||||||
|
|---|---|---|
|
||||||
|
| `test_model_manager.py` | ModelManager lifecycle (list, load, validate) | ✅ Well-tested |
|
||||||
|
| `test_feature_extractor.py` | All 4 feature extraction modes | ✅ Well-tested |
|
||||||
|
| `test_training_api.py` | Training API endpoints | ✅ Partial |
|
||||||
|
| `test_cloud_removal.py` | 7 cloud removal strategies | ✅ Well-tested |
|
||||||
|
| `test_cloud_training.py` | U-Net cloud removal training | ✅ Partial |
|
||||||
|
| `test_shapefile_api.py` | Shapefile overlay feature | ✅ Partial |
|
||||||
|
| `test_planetary_computer.py` | Planetary Computer STAC access | ✅ Well-tested |
|
||||||
|
| `test_new_features.py` | Recent feature releases | ✅ Partial |
|
||||||
|
| Jupyter Notebooks | Training & prediction workflows | ✅ Mix of unit/integration/notebooks |
|
||||||
|
|
||||||
|
**Coverage Notes**: Model management, feature extraction, and cloud removal well-tested; Dashboard UI, change detection, NDVI time-series mostly tested via notebooks.
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 📚 FILE REFERENCE MAP
|
||||||
|
|
||||||
|
### Core Execution
|
||||||
|
- `api_server.py` — Main FastAPI application (~4200 LOC)
|
||||||
|
- `train_module.py` — Training logic (data fetch → feature extraction → training)
|
||||||
|
- `run_prediction_new.py` — Prediction execution function
|
||||||
|
- `feature_extractor.py` — Unified feature extraction (4 modes)
|
||||||
|
- `model_manager.py` — Model lifecycle (load/save/validate/list)
|
||||||
|
- `cloud_removal.py` — Cloud removal strategies (7 methods)
|
||||||
|
- `report_generator.py` — HTML/PNG report auto-generation
|
||||||
|
- `generate_previews.py` — GeoTIFF → PNG conversion
|
||||||
|
|
||||||
|
### Data & Config
|
||||||
|
- `vietnam_provinces.py` — 32+ province lookup tables & bboxes
|
||||||
|
- `vietnam_provinces_merged.py` — Alternative 32-province variant
|
||||||
|
- `utils.py` — Geospatial utility functions
|
||||||
|
- `create_odc_metadata.py` — Legacy metadata generator
|
||||||
|
|
||||||
|
### Frontend
|
||||||
|
- `index.html` — Main dashboard hub (tab navigation)
|
||||||
|
- `training_interface.html` — Training configuration UI
|
||||||
|
- `prediction_interface.html` — Prediction configuration UI
|
||||||
|
- `batch_interface.html` — Batch processing (CSV upload)
|
||||||
|
- `ndvi_interface.html` — NDVI time-series visualization
|
||||||
|
- `dashboard.html` — Analytics & model performance dashboard
|
||||||
|
- `reports_interface.html` — Report management & viewing
|
||||||
|
- `change_detection_interface.html` — Change detection visualization
|
||||||
|
- `cloud_training_interface.html` — Cloud removal U-Net training
|
||||||
|
|
||||||
|
### Documentation
|
||||||
|
- `IMPLEMENTATION_SUMMARY.md` — Model manager & system overview
|
||||||
|
- `MODEL_MANAGER_GUIDE.md` — Complete model management guide
|
||||||
|
- `NDVI_FORECAST_METHODOLOGY.md` — NDVI algorithm documentation
|
||||||
|
- `CLOUD_TRAINING_GUIDE.md` — Cloud removal training guide
|
||||||
|
- `NDVI_PREDICTION_GUIDE.md` — NDVI prediction workflow
|
||||||
|
- `CLOUD_PROCESSING.md` — Cloud processing notes
|
||||||
|
- `UPDATE_SUMMARY.md` — Recent updates & features
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 🚀 BOOTSTRAP PROMPT FOR GEMINI
|
||||||
|
|
||||||
|
### System Context (Copy & Paste for Gemini)
|
||||||
|
|
||||||
|
```
|
||||||
|
You are assisting a remote-sensing land-classification project for Vietnam.
|
||||||
|
|
||||||
|
## ARCHITECTURE SNAPSHOT
|
||||||
|
- **Backend**: FastAPI (~4200 LOC, 70+ endpoints) for orchestrating training, prediction, batch, cloud removal, reporting
|
||||||
|
- **Data Source**: Microsoft Planetary Computer STAC API (Sentinel-2 L2A + Sentinel-1 RTC)
|
||||||
|
- **Training**: scikit-learn (XGBoost/RF/SVM/DT) + PyTorch (CNN/Swin-UNet/MobileNet)
|
||||||
|
- **Feature Extraction**: 4 modes (simple 3-feat / temporal 39-feat / extended 15-feat / odc 8-feat)
|
||||||
|
- **Cloud Removal**: 7 strategies (classic, temporal_only, median_composite, none, ML inpainting, deep U-Net)
|
||||||
|
- **Output**: GeoTIFF + PNG + HTML report + JSON metadata
|
||||||
|
|
||||||
|
## CORE FILES TO UNDERSTAND (Priority Order)
|
||||||
|
1. api_server.py — Main API server (training, prediction, batch, models, reports)
|
||||||
|
2. train_module.py — Training pipeline (data fetch → feature extraction → train → save)
|
||||||
|
3. feature_extractor.py — Unified feature extraction with auto mode detection
|
||||||
|
4. model_manager.py — Model lifecycle (load/save/validate/list)
|
||||||
|
5. cloud_removal.py — Cloud removal strategies (7 methods)
|
||||||
|
6. report_generator.py — Auto-generate HTML reports
|
||||||
|
7. run_prediction_new.py — Prediction execution
|
||||||
|
8. vietnam_provinces.py — Province lookup & bbox tables
|
||||||
|
|
||||||
|
## CRITICAL CONSTRAINTS & GOTCHAS
|
||||||
|
1. **Feature Mode Consistency**: Training & prediction MUST use same mode (simple/temporal/extended/odc)
|
||||||
|
→ Auto-detected from metadata JSON
|
||||||
|
→ Mismatch causes dimension error or accuracy degradation
|
||||||
|
|
||||||
|
2. **Planetary Computer Limits**:
|
||||||
|
→ Timeout if bbox >10km×10km OR time range >1 month OR max_scenes >12
|
||||||
|
→ Solution: subdivide bbox, reduce time window, limit scenes
|
||||||
|
|
||||||
|
3. **Cloud Strategy Selection**:
|
||||||
|
→ ≥5 scenes → use "median_composite" (best noise reduction)
|
||||||
|
→ 3-4 scenes → use "temporal_only" (fast temporal interp)
|
||||||
|
→ <3 scenes → use "none" (skip cloud removal)
|
||||||
|
|
||||||
|
4. **Radar Data Fallback**:
|
||||||
|
→ Sentinel-1 may be unavailable for some regions
|
||||||
|
→ System gracefully falls back to zeros (safe for all modes)
|
||||||
|
|
||||||
|
5. **Model Metadata**:
|
||||||
|
→ Always stored as `model_name_info.json` sidecar file
|
||||||
|
→ Contains: n_features, feature_mode, features list, accuracy, bbox, time_range
|
||||||
|
→ Missing metadata → system uses defaults (may be incorrect)
|
||||||
|
|
||||||
|
6. **Legacy Model (model_odc.joblib)**:
|
||||||
|
→ Hardcoded 39 temporal features
|
||||||
|
→ Metadata in model_odc_info.json
|
||||||
|
|
||||||
|
## REASONING CHECKLIST (before answering)
|
||||||
|
□ Is feature_mode consistent between train and prediction?
|
||||||
|
□ Is metadata.json present and correct?
|
||||||
|
□ Does bbox exceed 10km×10km? (Planetary Computer timeout risk)
|
||||||
|
□ Is cloud_removal_strategy appropriate for # of scenes?
|
||||||
|
□ Is Sentinel-1 data available for this region/date?
|
||||||
|
□ Is model a joblib (scikit-learn) or .pth (PyTorch) file?
|
||||||
|
□ Is GPU available for deep models (CNN, Swin-UNet)?
|
||||||
|
□ Does memory allow temporal feature extraction (39-feat)?
|
||||||
|
|
||||||
|
## RESPONSE FORMAT
|
||||||
|
- Always cite api_server.py endpoint, function name, or module being discussed
|
||||||
|
- Verify feature_mode & n_features from metadata JSON
|
||||||
|
- Suggest cloud_removal_strategy based on # of scenes available
|
||||||
|
- For unknown issues: offer alternative approaches (reduce bbox, cache results, use simpler model)
|
||||||
|
- Explain reasoning using checklist above
|
||||||
|
|
||||||
|
## DATA FLOW SUMMARY
|
||||||
|
Sentinel-2/S1 → [Cloud Remove] → [Feature Extract] → [Train/Predict] → [GeoTIFF + PNG + Report]
|
||||||
|
```
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 📞 QUICK REFERENCE CHECKLIST
|
||||||
|
|
||||||
|
### Before Troubleshooting Any Issue
|
||||||
|
- [ ] Check feature_mode consistency (metadata JSON)
|
||||||
|
- [ ] Verify metadata.json exists for the model
|
||||||
|
- [ ] Check Planetary Computer connectivity (`GET /api/network/check`)
|
||||||
|
- [ ] Review cloud_removal_method choice (≥5 scenes = median_composite)
|
||||||
|
- [ ] Confirm Sentinel-1 availability (or fallback to zeros if missing)
|
||||||
|
- [ ] Validate bbox size (≤10km×10km for safety)
|
||||||
|
- [ ] Check memory usage for temporal feature mode
|
||||||
|
- [ ] Verify GPU if using CNN/Swin-UNet models
|
||||||
|
|
||||||
|
### Common Issues & Solutions
|
||||||
|
| Issue | Likely Cause | Solution |
|
||||||
|
|---|---|---|
|
||||||
|
| Training timeout | Large bbox / long time / many scenes | Subdivide bbox, reduce time window, max_scenes ≤ 12 |
|
||||||
|
| Feature dimension mismatch | Different feature_mode between train & predict | Check metadata.json, ensure same mode |
|
||||||
|
| Low prediction accuracy | Cloud cover, poor training data, feature mode too simple | Use "temporal" mode, increase training data, try cloud removal |
|
||||||
|
| Out of memory | Temporal features + large region | Reduce resolution (20m), split into sub-tiles, increase RAM |
|
||||||
|
| Model not found | Wrong filename or model_train/ path issue | `GET /api/models/list` to verify, check file path |
|
||||||
|
| Planetary Computer error | Network issue or API rate limit | Check DNS, retry later, reduce concurrent requests |
|
||||||
|
| Cloud removal failing | Strategy not suitable for scene count | Try "none" or "median_composite" depending on scenes |
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
## 🎓 LEARNING RESOURCES IN REPO
|
||||||
|
|
||||||
|
- **Notebooks**: `01.train_ODC.ipynb`, `02.predict_ODC.ipynb`, `cloud_removal_train.ipynb`
|
||||||
|
- **Tests**: `test_*.py` files for unit test patterns
|
||||||
|
- **Docs**: All `*.md` files for detailed guides and methodology
|
||||||
|
- **Code Comments**: API server and modules heavily commented
|
||||||
|
|
||||||
|
---
|
||||||
|
|
||||||
|
**Generated**: March 26, 2026
|
||||||
|
**For Use By**: Gemini, Claude, GPT, or any AI system needing project context
|
||||||
|
**Maintainer**: Remote-Sensing Project Team
|
||||||
|
|
||||||
+296
-57
@@ -26,20 +26,20 @@ import socket
|
|||||||
import urllib.request
|
import urllib.request
|
||||||
|
|
||||||
# Import report generator
|
# Import report generator
|
||||||
from scripts.inference.report_generator import generate_training_report, generate_prediction_report
|
from report_generator import generate_training_report, generate_prediction_report
|
||||||
|
|
||||||
# Import Model Manager
|
# Import Model Manager
|
||||||
from core.model_manager import ModelManager, get_model_manager
|
from model_manager import ModelManager, get_model_manager
|
||||||
|
|
||||||
# Import Vietnam provinces data
|
# Import Vietnam provinces data
|
||||||
from core.vietnam_provinces import get_all_provinces, get_provinces_by_region, get_province_bbox, search_province
|
from vietnam_provinces import get_all_provinces, get_provinces_by_region, get_province_bbox, search_province
|
||||||
from core.vietnam_provinces_merged import (
|
from vietnam_provinces_merged import (
|
||||||
get_all_provinces_32, get_provinces_by_region_32, get_province_bbox_32,
|
get_all_provinces_32, get_provinces_by_region_32, get_province_bbox_32,
|
||||||
search_province_32, get_merged_info, get_provinces_statistics
|
search_province_32, get_merged_info, get_provinces_statistics
|
||||||
)
|
)
|
||||||
|
|
||||||
# Import cloud removal module
|
# Import cloud removal module
|
||||||
from core.cloud_removal import process_cloud_removal, get_available_methods
|
from cloud_removal import process_cloud_removal, get_available_methods
|
||||||
|
|
||||||
# Import planetary computer libraries (conditional)
|
# Import planetary computer libraries (conditional)
|
||||||
try:
|
try:
|
||||||
@@ -225,6 +225,9 @@ class PredictionConfig(BaseModel):
|
|||||||
cloud_removal_method: str = "classic"
|
cloud_removal_method: str = "classic"
|
||||||
cloud_removal_model: Optional[str] = None # Optional: .pth model filename for deep learning cloud removal
|
cloud_removal_model: Optional[str] = None # Optional: .pth model filename for deep learning cloud removal
|
||||||
|
|
||||||
|
# Shapefile overlay for visualization
|
||||||
|
shapefile_overlay: Optional[str] = None # Path to shapefile for overlaying boundaries
|
||||||
|
|
||||||
|
|
||||||
class TrainingStatus(BaseModel):
|
class TrainingStatus(BaseModel):
|
||||||
"""Trạng thái training"""
|
"""Trạng thái training"""
|
||||||
@@ -286,6 +289,7 @@ class PredictionWithNDVIConfig(BaseModel):
|
|||||||
export_classification: bool = True # Export classification raster
|
export_classification: bool = True # Export classification raster
|
||||||
cloud_removal_method: str = "classic"
|
cloud_removal_method: str = "classic"
|
||||||
cloud_removal_model: Optional[str] = None # Optional: .pth model filename for deep learning cloud removal
|
cloud_removal_model: Optional[str] = None # Optional: .pth model filename for deep learning cloud removal
|
||||||
|
shapefile_overlay: Optional[str] = None # Path to shapefile for overlaying boundaries
|
||||||
|
|
||||||
|
|
||||||
class CloudRemovalTrainingConfig(BaseModel):
|
class CloudRemovalTrainingConfig(BaseModel):
|
||||||
@@ -507,57 +511,6 @@ async def get_cloud_removal_methods():
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
@app.get("/api/ndvi-forecast/models")
|
|
||||||
async def list_ndvi_forecast_models():
|
|
||||||
"""Liệt kê các NDVI forecast models đã train"""
|
|
||||||
model_dir = Path("ndvi_forecast_model")
|
|
||||||
if not model_dir.exists():
|
|
||||||
return {"models": [], "count": 0}
|
|
||||||
|
|
||||||
models = []
|
|
||||||
# Search for all models
|
|
||||||
for model_file in list(model_dir.rglob("*.pth")) + list(model_dir.rglob("*.joblib")):
|
|
||||||
try:
|
|
||||||
import json
|
|
||||||
|
|
||||||
# Try to load metadata from .json sidecar file first
|
|
||||||
metadata_file = model_file.with_name(model_file.stem + "_info.json")
|
|
||||||
if metadata_file.exists():
|
|
||||||
try:
|
|
||||||
with open(metadata_file, 'r') as f:
|
|
||||||
metadata = json.load(f)
|
|
||||||
|
|
||||||
models.append({
|
|
||||||
"filename": model_file.name,
|
|
||||||
"path": str(model_file),
|
|
||||||
"model_type": metadata.get('model_type', 'Unknown'),
|
|
||||||
"target": metadata.get('target', 'NDVI'),
|
|
||||||
"rmse": metadata.get('rmse', 0),
|
|
||||||
"mae": metadata.get('mae', 0),
|
|
||||||
"epoch": metadata.get('epoch', 0),
|
|
||||||
"created": model_file.stat().st_mtime,
|
|
||||||
"size_mb": model_file.stat().st_size / (1024 * 1024),
|
|
||||||
})
|
|
||||||
continue
|
|
||||||
except Exception as e:
|
|
||||||
print(f"Error reading JSON {metadata_file}: {e}")
|
|
||||||
|
|
||||||
# Fallback for models without metadata
|
|
||||||
models.append({
|
|
||||||
"filename": model_file.name,
|
|
||||||
"path": str(model_file),
|
|
||||||
"model_type": "Unknown",
|
|
||||||
"created": model_file.stat().st_mtime,
|
|
||||||
"size_mb": model_file.stat().st_size / (1024 * 1024)
|
|
||||||
})
|
|
||||||
except Exception as e:
|
|
||||||
print(f"Error loading model info for {model_file}: {e}")
|
|
||||||
|
|
||||||
# Sort by creation time (newest first)
|
|
||||||
models.sort(key=lambda x: x['created'], reverse=True)
|
|
||||||
return {"models": models, "count": len(models)}
|
|
||||||
|
|
||||||
|
|
||||||
@app.get("/api/cloud-removal/models")
|
@app.get("/api/cloud-removal/models")
|
||||||
async def list_cloud_removal_models():
|
async def list_cloud_removal_models():
|
||||||
"""Liệt kê các cloud removal models đã train"""
|
"""Liệt kê các cloud removal models đã train"""
|
||||||
@@ -1279,6 +1232,87 @@ async def list_training_files():
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
@app.get("/api/overlay/shapefiles")
|
||||||
|
async def list_overlay_shapefiles():
|
||||||
|
"""Liệt kê các shapefile có sẵn cho overlay trên prediction"""
|
||||||
|
overlay_dirs = ["region", "ChauThanh", "ThuanHoa"]
|
||||||
|
shapefiles = []
|
||||||
|
|
||||||
|
for overlay_dir in overlay_dirs:
|
||||||
|
dir_path = Path(overlay_dir)
|
||||||
|
if not dir_path.exists():
|
||||||
|
continue
|
||||||
|
|
||||||
|
# Find all .shp files in this directory and subdirectories
|
||||||
|
for shp_file in dir_path.rglob("*.shp"):
|
||||||
|
try:
|
||||||
|
file_size = shp_file.stat().st_size
|
||||||
|
file_modified = datetime.fromtimestamp(shp_file.stat().st_mtime).isoformat()
|
||||||
|
|
||||||
|
# Try to read shapefile to get feature count and bbox
|
||||||
|
import geopandas as gpd
|
||||||
|
gdf = gpd.read_file(str(shp_file))
|
||||||
|
feature_count = len(gdf)
|
||||||
|
|
||||||
|
# Calculate bbox (always in EPSG:4326 for consistency)
|
||||||
|
bbox = None
|
||||||
|
if not gdf.empty and gdf.crs:
|
||||||
|
try:
|
||||||
|
# Reproject to EPSG:4326 if needed
|
||||||
|
if gdf.crs != "EPSG:4326":
|
||||||
|
gdf_4326 = gdf.to_crs("EPSG:4326")
|
||||||
|
else:
|
||||||
|
gdf_4326 = gdf
|
||||||
|
|
||||||
|
# Get total bounds [minx, miny, maxx, maxy]
|
||||||
|
bounds = gdf_4326.total_bounds
|
||||||
|
if len(bounds) == 4:
|
||||||
|
bbox = [
|
||||||
|
float(bounds[0]), # min_lon
|
||||||
|
float(bounds[1]), # min_lat
|
||||||
|
float(bounds[2]), # max_lon
|
||||||
|
float(bounds[3]) # max_lat
|
||||||
|
]
|
||||||
|
except Exception as bbox_error:
|
||||||
|
print(f"[WARNING] Cannot calculate bbox for {shp_file}: {bbox_error}")
|
||||||
|
|
||||||
|
# Get relative path from workspace root
|
||||||
|
relative_path = str(shp_file)
|
||||||
|
|
||||||
|
shapefiles.append({
|
||||||
|
"filename": shp_file.name,
|
||||||
|
"path": relative_path,
|
||||||
|
"directory": overlay_dir,
|
||||||
|
"size_bytes": file_size,
|
||||||
|
"size_mb": round(file_size / 1024 / 1024, 2),
|
||||||
|
"modified": file_modified,
|
||||||
|
"feature_count": feature_count,
|
||||||
|
"crs": str(gdf.crs) if gdf.crs else "Unknown",
|
||||||
|
"bbox": bbox # [min_lon, min_lat, max_lon, max_lat] in EPSG:4326
|
||||||
|
})
|
||||||
|
except Exception as e:
|
||||||
|
# If cannot read shapefile, just add basic info
|
||||||
|
file_size = shp_file.stat().st_size
|
||||||
|
file_modified = datetime.fromtimestamp(shp_file.stat().st_mtime).isoformat()
|
||||||
|
relative_path = str(shp_file)
|
||||||
|
|
||||||
|
shapefiles.append({
|
||||||
|
"filename": shp_file.name,
|
||||||
|
"path": relative_path,
|
||||||
|
"directory": overlay_dir,
|
||||||
|
"size_bytes": file_size,
|
||||||
|
"size_mb": round(file_size / 1024 / 1024, 2),
|
||||||
|
"modified": file_modified,
|
||||||
|
"error": f"Cannot read shapefile: {str(e)}"
|
||||||
|
})
|
||||||
|
|
||||||
|
return {
|
||||||
|
"shapefiles": shapefiles,
|
||||||
|
"count": len(shapefiles),
|
||||||
|
"directories": overlay_dirs
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
@app.get("/api/training/shapefile/{filename}/labels")
|
@app.get("/api/training/shapefile/{filename}/labels")
|
||||||
async def get_shapefile_labels(filename: str):
|
async def get_shapefile_labels(filename: str):
|
||||||
"""Lấy các label từ một shapefile cụ thể"""
|
"""Lấy các label từ một shapefile cụ thể"""
|
||||||
@@ -1780,6 +1814,85 @@ def update_progress(message: str):
|
|||||||
print(f"[PROGRESS] {message}")
|
print(f"[PROGRESS] {message}")
|
||||||
|
|
||||||
|
|
||||||
|
def rasterize_shapefile_overlay(shapefile_path, reference_raster, boundary_value=255):
|
||||||
|
"""
|
||||||
|
Rasterize shapefile boundaries to overlay on prediction result.
|
||||||
|
|
||||||
|
Args:
|
||||||
|
shapefile_path: Path to shapefile
|
||||||
|
reference_raster: xarray DataArray to match dimensions and CRS
|
||||||
|
boundary_value: Value to use for boundaries (default 255 for white)
|
||||||
|
|
||||||
|
Returns:
|
||||||
|
numpy array with boundaries, same shape as reference_raster
|
||||||
|
"""
|
||||||
|
try:
|
||||||
|
import geopandas as gpd
|
||||||
|
from rasterio.features import rasterize
|
||||||
|
import numpy as np
|
||||||
|
|
||||||
|
# Read shapefile
|
||||||
|
gdf = gpd.read_file(shapefile_path)
|
||||||
|
print(f"[OVERLAY] Loaded shapefile with {len(gdf)} features, CRS: {gdf.crs}")
|
||||||
|
|
||||||
|
# Ensure CRS matches
|
||||||
|
target_crs = reference_raster.rio.crs
|
||||||
|
if gdf.crs != target_crs:
|
||||||
|
print(f"[OVERLAY] Reprojecting from {gdf.crs} to {target_crs}")
|
||||||
|
gdf = gdf.to_crs(target_crs)
|
||||||
|
|
||||||
|
# Get raster dimensions and transform
|
||||||
|
height, width = reference_raster.shape
|
||||||
|
transform = reference_raster.rio.transform()
|
||||||
|
print(f"[OVERLAY] Raster dimensions: {height}x{width}")
|
||||||
|
print(f"[OVERLAY] Transform: {transform}")
|
||||||
|
|
||||||
|
# Calculate appropriate buffer size based on pixel resolution
|
||||||
|
# Get pixel size from transform (transform[0] is x resolution)
|
||||||
|
pixel_size = abs(transform[0]) # in CRS units
|
||||||
|
# Very thin boundary - only 0.2 pixels wide for 1px line
|
||||||
|
buffer_distance = pixel_size * 0.2
|
||||||
|
print(f"[OVERLAY] Pixel size: {pixel_size}, Buffer distance: {buffer_distance} (thin 1px line)")
|
||||||
|
|
||||||
|
# Create boundary geometries with minimal buffering
|
||||||
|
boundary_geoms = []
|
||||||
|
for idx, geom in enumerate(gdf.geometry):
|
||||||
|
if geom is not None and geom.is_valid:
|
||||||
|
# Get boundary of each polygon
|
||||||
|
boundary = geom.boundary
|
||||||
|
if boundary is not None:
|
||||||
|
# Minimal buffer for 1-pixel thin line
|
||||||
|
buffered = boundary.buffer(buffer_distance)
|
||||||
|
boundary_geoms.append((buffered, boundary_value))
|
||||||
|
|
||||||
|
if not boundary_geoms:
|
||||||
|
print(f"[WARNING] No valid boundary geometries found in {shapefile_path}")
|
||||||
|
return np.zeros((height, width), dtype=np.uint8)
|
||||||
|
|
||||||
|
print(f"[OVERLAY] Rasterizing {len(boundary_geoms)} boundaries...")
|
||||||
|
|
||||||
|
# Rasterize boundaries
|
||||||
|
boundary_mask = rasterize(
|
||||||
|
shapes=boundary_geoms,
|
||||||
|
out_shape=(height, width),
|
||||||
|
transform=transform,
|
||||||
|
fill=0, # Background
|
||||||
|
dtype=np.uint8
|
||||||
|
)
|
||||||
|
|
||||||
|
boundary_count = np.count_nonzero(boundary_mask)
|
||||||
|
print(f"[OVERLAY] Boundary pixels: {boundary_count} / {height*width} ({boundary_count/(height*width)*100:.2f}%)")
|
||||||
|
|
||||||
|
if boundary_count == 0:
|
||||||
|
print(f"[OVERLAY WARNING] No boundary pixels were rasterized! Check CRS and geometry overlap.")
|
||||||
|
|
||||||
|
return boundary_mask
|
||||||
|
|
||||||
|
except Exception as e:
|
||||||
|
print(f"[ERROR] Failed to rasterize shapefile {shapefile_path}: {e}")
|
||||||
|
return None
|
||||||
|
|
||||||
|
|
||||||
def update_prediction_progress(message: str):
|
def update_prediction_progress(message: str):
|
||||||
"""Cập nhật prediction progress message"""
|
"""Cập nhật prediction progress message"""
|
||||||
global prediction_status
|
global prediction_status
|
||||||
@@ -2134,6 +2247,33 @@ async def run_prediction(config: PredictionConfig):
|
|||||||
|
|
||||||
prediction_da.rio.to_raster(str(output_file), driver="GTiff")
|
prediction_da.rio.to_raster(str(output_file), driver="GTiff")
|
||||||
|
|
||||||
|
# ============ SHAPEFILE OVERLAY ============
|
||||||
|
overlay_mask = None
|
||||||
|
if config.shapefile_overlay:
|
||||||
|
prediction_status["progress"] = f"Đang overlay shapefile: {config.shapefile_overlay}..."
|
||||||
|
print(f"[OVERLAY] Shapefile overlay requested: {config.shapefile_overlay}")
|
||||||
|
|
||||||
|
# Validate shapefile path exists
|
||||||
|
shapefile_path = Path(config.shapefile_overlay)
|
||||||
|
if not shapefile_path.exists():
|
||||||
|
print(f"[OVERLAY ERROR] Shapefile not found: {shapefile_path}")
|
||||||
|
print(f"[OVERLAY ERROR] Absolute path: {shapefile_path.absolute()}")
|
||||||
|
print(f"[OVERLAY ERROR] Current working directory: {Path.cwd()}")
|
||||||
|
else:
|
||||||
|
print(f"[OVERLAY] Shapefile exists: {shapefile_path.absolute()}")
|
||||||
|
try:
|
||||||
|
overlay_mask = rasterize_shapefile_overlay(str(shapefile_path), prediction_da)
|
||||||
|
if overlay_mask is not None and np.count_nonzero(overlay_mask) > 0:
|
||||||
|
print(f"[OVERLAY] Successfully rasterized shapefile boundaries ({np.count_nonzero(overlay_mask)} pixels)")
|
||||||
|
else:
|
||||||
|
print(f"[OVERLAY WARNING] Shapefile rasterized but no boundary pixels found")
|
||||||
|
except Exception as overlay_error:
|
||||||
|
print(f"[OVERLAY ERROR] Exception: {overlay_error}")
|
||||||
|
import traceback
|
||||||
|
traceback.print_exc()
|
||||||
|
else:
|
||||||
|
print(f"[OVERLAY] No shapefile overlay requested")
|
||||||
|
|
||||||
# Generate PNG preview for web display
|
# Generate PNG preview for web display
|
||||||
prediction_status["progress"] = "Đang tạo PNG preview..."
|
prediction_status["progress"] = "Đang tạo PNG preview..."
|
||||||
png_file = output_dir / f"prediction_{timestamp}.png"
|
png_file = output_dir / f"prediction_{timestamp}.png"
|
||||||
@@ -2148,6 +2288,42 @@ async def run_prediction(config: PredictionConfig):
|
|||||||
|
|
||||||
# Plot prediction with colormap
|
# Plot prediction with colormap
|
||||||
im = ax.imshow(predictions_2d, cmap='tab20', interpolation='nearest')
|
im = ax.imshow(predictions_2d, cmap='tab20', interpolation='nearest')
|
||||||
|
|
||||||
|
# Overlay shapefile boundaries if available
|
||||||
|
if overlay_mask is not None and np.count_nonzero(overlay_mask) > 0:
|
||||||
|
print(f"[PNG OVERLAY] Overlaying {np.count_nonzero(overlay_mask)} boundary pixels")
|
||||||
|
|
||||||
|
# Create a mask for boundaries (where overlay_mask > 0)
|
||||||
|
boundary_mask = overlay_mask > 0
|
||||||
|
|
||||||
|
# Method: Direct overlay with high-contrast colors
|
||||||
|
# Create RGBA overlay image
|
||||||
|
overlay_rgba = np.zeros((*predictions_2d.shape, 4))
|
||||||
|
overlay_rgba[boundary_mask, 0] = 1.0 # Red = 1.0 (white)
|
||||||
|
overlay_rgba[boundary_mask, 1] = 1.0 # Green = 1.0 (white)
|
||||||
|
overlay_rgba[boundary_mask, 2] = 1.0 # Blue = 1.0 (white)
|
||||||
|
overlay_rgba[boundary_mask, 3] = 1.0 # Alpha = 1.0 (fully opaque)
|
||||||
|
|
||||||
|
# Overlay on top of prediction
|
||||||
|
ax.imshow(overlay_rgba, interpolation='nearest')
|
||||||
|
|
||||||
|
# Also add a black outline for better contrast
|
||||||
|
from scipy import ndimage
|
||||||
|
boundary_dilated = ndimage.binary_dilation(boundary_mask, iterations=1)
|
||||||
|
boundary_outline = boundary_dilated & ~boundary_mask
|
||||||
|
|
||||||
|
outline_rgba = np.zeros((*predictions_2d.shape, 4))
|
||||||
|
outline_rgba[boundary_outline, 0] = 0.0 # Black outline
|
||||||
|
outline_rgba[boundary_outline, 1] = 0.0
|
||||||
|
outline_rgba[boundary_outline, 2] = 0.0
|
||||||
|
outline_rgba[boundary_outline, 3] = 0.8
|
||||||
|
|
||||||
|
ax.imshow(outline_rgba, interpolation='nearest')
|
||||||
|
|
||||||
|
print(f"[PNG OVERLAY] Added shapefile boundaries to visualization (direct overlay method)")
|
||||||
|
else:
|
||||||
|
print(f"[PNG OVERLAY] No overlay mask or empty mask (pixels: {np.count_nonzero(overlay_mask) if overlay_mask is not None else 0})")
|
||||||
|
|
||||||
ax.set_title(f'Land Classification - {timestamp}', fontsize=16, fontweight='bold', pad=20)
|
ax.set_title(f'Land Classification - {timestamp}', fontsize=16, fontweight='bold', pad=20)
|
||||||
ax.set_xlabel('X (pixels)', fontsize=11)
|
ax.set_xlabel('X (pixels)', fontsize=11)
|
||||||
ax.set_ylabel('Y (pixels)', fontsize=11)
|
ax.set_ylabel('Y (pixels)', fontsize=11)
|
||||||
@@ -2227,7 +2403,9 @@ async def run_prediction(config: PredictionConfig):
|
|||||||
"n_features": features.shape[1],
|
"n_features": features.shape[1],
|
||||||
"feature_mode": feature_mode,
|
"feature_mode": feature_mode,
|
||||||
"used_radar": use_radar,
|
"used_radar": use_radar,
|
||||||
"model_used": config.model_filename
|
"model_used": config.model_filename,
|
||||||
|
"shapefile_overlay": config.shapefile_overlay,
|
||||||
|
"overlay_applied": overlay_mask is not None
|
||||||
}
|
}
|
||||||
|
|
||||||
# Auto generate prediction report
|
# Auto generate prediction report
|
||||||
@@ -4406,8 +4584,69 @@ async def predict_with_ndvi(config: PredictionWithNDVIConfig, background_tasks:
|
|||||||
import matplotlib.pyplot as plt
|
import matplotlib.pyplot as plt
|
||||||
from matplotlib.patches import Patch
|
from matplotlib.patches import Patch
|
||||||
|
|
||||||
|
# ============ SHAPEFILE OVERLAY ============
|
||||||
|
overlay_mask = None
|
||||||
|
if config.shapefile_overlay:
|
||||||
|
print(f"[OVERLAY] Shapefile overlay requested: {config.shapefile_overlay}")
|
||||||
|
|
||||||
|
# Validate shapefile path exists
|
||||||
|
shapefile_path = Path(config.shapefile_overlay)
|
||||||
|
if not shapefile_path.exists():
|
||||||
|
print(f"[OVERLAY ERROR] Shapefile not found: {shapefile_path}")
|
||||||
|
print(f"[OVERLAY ERROR] Absolute path: {shapefile_path.absolute()}")
|
||||||
|
else:
|
||||||
|
print(f"[OVERLAY] Shapefile exists: {shapefile_path.absolute()}")
|
||||||
|
try:
|
||||||
|
# Import rioxarray for rio accessor
|
||||||
|
import rioxarray
|
||||||
|
|
||||||
|
# Create temporary DataArray for rasterization
|
||||||
|
temp_da = xr.DataArray(
|
||||||
|
prediction_raster,
|
||||||
|
coords={
|
||||||
|
"y": np.linspace(bbox[3], bbox[1], height),
|
||||||
|
"x": np.linspace(bbox[0], bbox[2], width)
|
||||||
|
},
|
||||||
|
dims=["y", "x"]
|
||||||
|
)
|
||||||
|
temp_da.rio.write_crs("EPSG:4326", inplace=True)
|
||||||
|
temp_da.rio.write_transform(transform, inplace=True)
|
||||||
|
|
||||||
|
overlay_mask = rasterize_shapefile_overlay(str(shapefile_path), temp_da)
|
||||||
|
if overlay_mask is not None and np.count_nonzero(overlay_mask) > 0:
|
||||||
|
print(f"[OVERLAY] Successfully rasterized shapefile boundaries ({np.count_nonzero(overlay_mask)} pixels)")
|
||||||
|
else:
|
||||||
|
print(f"[OVERLAY WARNING] Shapefile rasterized but no boundary pixels found")
|
||||||
|
except Exception as overlay_error:
|
||||||
|
print(f"[OVERLAY ERROR] Exception: {overlay_error}")
|
||||||
|
import traceback
|
||||||
|
traceback.print_exc()
|
||||||
|
else:
|
||||||
|
print(f"[OVERLAY] No shapefile overlay requested")
|
||||||
|
|
||||||
fig, ax = plt.subplots(figsize=(14, 10), dpi=150)
|
fig, ax = plt.subplots(figsize=(14, 10), dpi=150)
|
||||||
im = ax.imshow(prediction_raster, cmap='tab20', interpolation='nearest')
|
im = ax.imshow(prediction_raster, cmap='tab20', interpolation='nearest')
|
||||||
|
|
||||||
|
# Overlay shapefile boundaries if available
|
||||||
|
if overlay_mask is not None and np.count_nonzero(overlay_mask) > 0:
|
||||||
|
print(f"[PNG OVERLAY] Overlaying {np.count_nonzero(overlay_mask)} boundary pixels (1px thin line)")
|
||||||
|
|
||||||
|
# Create a mask for boundaries
|
||||||
|
boundary_mask = overlay_mask > 0
|
||||||
|
|
||||||
|
# Create RGBA overlay image - thin 1px white line only
|
||||||
|
overlay_rgba = np.zeros((*prediction_raster.shape, 4))
|
||||||
|
overlay_rgba[boundary_mask, 0] = 1.0 # White (R=1)
|
||||||
|
overlay_rgba[boundary_mask, 1] = 1.0 # White (G=1)
|
||||||
|
overlay_rgba[boundary_mask, 2] = 1.0 # White (B=1)
|
||||||
|
overlay_rgba[boundary_mask, 3] = 1.0 # Fully opaque
|
||||||
|
|
||||||
|
ax.imshow(overlay_rgba, interpolation='nearest')
|
||||||
|
|
||||||
|
print(f"[PNG OVERLAY] Added thin 1px shapefile boundaries to visualization")
|
||||||
|
else:
|
||||||
|
print(f"[PNG OVERLAY] No overlay mask or empty mask")
|
||||||
|
|
||||||
ax.set_title(f'Land Classification - {timestamp}', fontsize=16, fontweight='bold', pad=20)
|
ax.set_title(f'Land Classification - {timestamp}', fontsize=16, fontweight='bold', pad=20)
|
||||||
ax.set_xlabel('X (pixels)', fontsize=11)
|
ax.set_xlabel('X (pixels)', fontsize=11)
|
||||||
ax.set_ylabel('Y (pixels)', fontsize=11)
|
ax.set_ylabel('Y (pixels)', fontsize=11)
|
||||||
|
|||||||
@@ -0,0 +1,354 @@
|
|||||||
|
import matplotlib.pyplot as plt
|
||||||
|
|
||||||
|
# Common imports and settings
|
||||||
|
import os, sys
|
||||||
|
os.environ['USE_PYGEOS'] = '0'
|
||||||
|
from IPython.display import Markdown
|
||||||
|
import pandas as pd
|
||||||
|
pd.set_option("display.max_rows", None)
|
||||||
|
import xarray as xr
|
||||||
|
|
||||||
|
# Datacube
|
||||||
|
import datacube
|
||||||
|
from datacube.utils.rio import configure_s3_access
|
||||||
|
from datacube.utils import masking
|
||||||
|
from datacube.utils.cog import write_cog
|
||||||
|
# https://github.com/GeoscienceAustralia/dea-notebooks/tree/develop/Tools
|
||||||
|
from dea_tools.plotting import display_map, rgb
|
||||||
|
from dea_tools.datahandling import mostcommon_crs
|
||||||
|
|
||||||
|
# EASI defaults
|
||||||
|
easinotebooksrepo = '/home/jovyan/easi-notebooks'
|
||||||
|
if easinotebooksrepo not in sys.path: sys.path.append(easinotebooksrepo)
|
||||||
|
from easi_tools import EasiDefaults, xarray_object_size, notebook_utils, unset_cachingproxy
|
||||||
|
from easi_tools.load_s2l2a import load_s2l2a_with_offset
|
||||||
|
from dask.distributed import progress
|
||||||
|
|
||||||
|
# Data tools
|
||||||
|
import numpy as np
|
||||||
|
from datetime import datetime
|
||||||
|
|
||||||
|
# Datacube
|
||||||
|
from datacube.utils import masking # https://github.com/opendatacube/datacube-core/blob/develop/datacube/utils/masking.py
|
||||||
|
from odc.algo import enum_to_bool # https://github.com/opendatacube/odc-algo/blob/main/odc/algo/_masking.py
|
||||||
|
from odc.algo import xr_reproject # https://github.com/opendatacube/odc-algo/blob/main/odc/algo/_warp.py
|
||||||
|
from datacube.utils.geometry import GeoBox, box # https://github.com/opendatacube/datacube-core/blob/develop/datacube/utils/geometry/_base.py
|
||||||
|
|
||||||
|
# Holoviews, Datashader and Bokeh
|
||||||
|
import hvplot.pandas
|
||||||
|
import hvplot.xarray
|
||||||
|
import holoviews as hv
|
||||||
|
import panel as pn
|
||||||
|
import colorcet as cc
|
||||||
|
import cartopy.crs as ccrs
|
||||||
|
from datashader import reductions
|
||||||
|
from holoviews import opts
|
||||||
|
from utils import load_data_geo
|
||||||
|
import rasterio
|
||||||
|
import rioxarray
|
||||||
|
# import geoviews as gv
|
||||||
|
# from holoviews.operation.datashader import rasterize
|
||||||
|
hv.extension('bokeh', logo=False)
|
||||||
|
|
||||||
|
from deafrica_tools.bandindices import calculate_indices
|
||||||
|
from sklearn.ensemble import RandomForestClassifier
|
||||||
|
from sklearn.model_selection import train_test_split
|
||||||
|
from sklearn.metrics import accuracy_score, classification_report
|
||||||
|
from sklearn.preprocessing import LabelEncoder
|
||||||
|
|
||||||
|
from sklearn.pipeline import Pipeline
|
||||||
|
from sklearn.ensemble import RandomForestClassifier
|
||||||
|
from sklearn.impute import SimpleImputer
|
||||||
|
from sklearn.preprocessing import StandardScaler
|
||||||
|
from sklearn.model_selection import GridSearchCV
|
||||||
|
from sklearn.model_selection import train_test_split
|
||||||
|
from sklearn.metrics import accuracy_score
|
||||||
|
from shapely.geometry import Point, Polygon
|
||||||
|
import geopandas as gpd
|
||||||
|
from pyproj import CRS
|
||||||
|
from matplotlib.colors import ListedColormap
|
||||||
|
from holoviews import opts
|
||||||
|
from datashader import reductions
|
||||||
|
from bokeh.models.tickers import FixedTicker
|
||||||
|
from rioxarray.merge import merge_arrays
|
||||||
|
|
||||||
|
import joblib
|
||||||
|
|
||||||
|
|
||||||
|
def load_data(dc, date_range, longtitude_range, latitude_range):
|
||||||
|
product = 's2_l2a'
|
||||||
|
query = {
|
||||||
|
'product': product, # Product name
|
||||||
|
'x': longtitude_range, # "x" axis bounds
|
||||||
|
'y': latitude_range, # "y" axis bounds
|
||||||
|
'time': date_range, # Any parsable date strings
|
||||||
|
}
|
||||||
|
native_crs = notebook_utils.mostcommon_crs(dc, query)
|
||||||
|
print(f'Most common native CRS: {native_crs}')
|
||||||
|
measurements = [ 'red', 'nir', 'scl']
|
||||||
|
|
||||||
|
load_params = {
|
||||||
|
'measurements': measurements, # Selected measurement or alias names
|
||||||
|
'output_crs': native_crs, # Target EPSG code
|
||||||
|
'resolution': (-10, 10), # Target resolution
|
||||||
|
'group_by': 'solar_day', # Scene grouping
|
||||||
|
'dask_chunks': {'x': 2048, 'y': 2048}, # Dask chunks
|
||||||
|
}
|
||||||
|
data = load_s2l2a_with_offset(
|
||||||
|
dc,
|
||||||
|
query | load_params # Combine the two dicts that contain our search and load parameters
|
||||||
|
)
|
||||||
|
return data
|
||||||
|
|
||||||
|
|
||||||
|
def mask_clean(data):
|
||||||
|
flag_name = 'scl'
|
||||||
|
flag_desc = masking.describe_variable_flags(data[flag_name]) # Pandas dataframe
|
||||||
|
display(flag_desc)
|
||||||
|
display(flag_desc.loc['qa'].values[1])
|
||||||
|
# Create a "data quality" Mask layer
|
||||||
|
flags_def = flag_desc.loc['qa'].values[1]
|
||||||
|
good_pixel_flags = [flags_def[str(i)] for i in [2, 4, 5, 6]] # To pass strings to enum_to_bool()
|
||||||
|
|
||||||
|
# enum_to_bool calculates the pixel-wise "or" of each set of pixels given by good_pixel_flags
|
||||||
|
# 1 = good data
|
||||||
|
# 0 = "bad" data
|
||||||
|
good_pixel_mask = enum_to_bool(data[flag_name], good_pixel_flags)
|
||||||
|
data_layer_names = [x for x in data.data_vars if x != 'scl']
|
||||||
|
# Apply good pixel mask to blue, green, red and nir.
|
||||||
|
result = data[data_layer_names].where(good_pixel_mask).persist()
|
||||||
|
return result
|
||||||
|
|
||||||
|
|
||||||
|
def fill_nan(ndvi, time_split):
|
||||||
|
rs = []
|
||||||
|
for times in time_split:
|
||||||
|
tmp = ndvi.sel(time=times)
|
||||||
|
fill_ds = tmp.sel(time=times).bfill(dim='time')
|
||||||
|
fill_ds = fill_ds.sel(time=times).ffill(dim='time')
|
||||||
|
rs.append(fill_ds)
|
||||||
|
merged_ndvi = xr.concat([i for i in rs], dim="time")
|
||||||
|
fill_m = merged_ndvi.bfill(dim="time")
|
||||||
|
fill_m = fill_m.ffill(dim="time")
|
||||||
|
return fill_m
|
||||||
|
|
||||||
|
|
||||||
|
def load_train_data(train_path):
|
||||||
|
train = load_data_geo(train_path)
|
||||||
|
return train
|
||||||
|
|
||||||
|
|
||||||
|
def load_sen1(name_vh, name_vv):
|
||||||
|
dsvv = rioxarray.open_rasterio(name_vv)
|
||||||
|
dsvh = rioxarray.open_rasterio(name_vh)
|
||||||
|
return dsvh, dsvv
|
||||||
|
|
||||||
|
|
||||||
|
def get_data_sen1_and_sen2(train, average_ndvi, dsvh, dsvv):
|
||||||
|
loaded_datasets = {}
|
||||||
|
for idx, point in train.iterrows():
|
||||||
|
key = f"point_{idx + 1}"
|
||||||
|
try:
|
||||||
|
ndvi_data = average_ndvi.sel(x=point.geometry.x, y=point.geometry.y, method='nearest').values
|
||||||
|
vh_data = dsvh.sel(x=point.geometry.x, y=point.geometry.y, method='nearest').values
|
||||||
|
vv_data = dsvv.sel(x=point.geometry.x, y=point.geometry.y, method='nearest').values
|
||||||
|
loaded_datasets[key] = {
|
||||||
|
"data": np.concatenate((ndvi_data, vh_data, vv_data)),
|
||||||
|
"label": point.HT_code
|
||||||
|
}
|
||||||
|
except Exception as e:
|
||||||
|
# loaded_datasets[key] = None
|
||||||
|
print(e)
|
||||||
|
return loaded_datasets
|
||||||
|
|
||||||
|
|
||||||
|
def split_train_data(train, label_mapping, datasets):
|
||||||
|
label_encoder = LabelEncoder()
|
||||||
|
|
||||||
|
# Fit and transform the labels
|
||||||
|
labels = train.Hientrang.values
|
||||||
|
numeric_labels = label_encoder.fit_transform([label_mapping[label] for label in labels])
|
||||||
|
X = []
|
||||||
|
x_new = []
|
||||||
|
lb_new = []
|
||||||
|
for k, v in datasets.items():
|
||||||
|
X.append(v)
|
||||||
|
for i in range(len(X)):
|
||||||
|
if X[i] is not None:
|
||||||
|
x_new.append(X[i]["data"])
|
||||||
|
lb_new.append(numeric_labels[i])
|
||||||
|
X_train, X_temp, y_train, y_temp= train_test_split(x_new, lb_new, test_size=0.4, random_state=42)
|
||||||
|
X_val, X_test, y_val, y_test = train_test_split(X_temp, y_temp, test_size=0.5, random_state=42)
|
||||||
|
return X_train, X_val, X_test, y_train, y_val, y_test
|
||||||
|
|
||||||
|
|
||||||
|
def train_with_rf(X_train, X_val, y_train, y_val):
|
||||||
|
# Takes 1-2 minutes to complete
|
||||||
|
|
||||||
|
# Tạo RandomForestClassifier mặc định để sử dụng làm mô hình ban đầu trong pipeline
|
||||||
|
base_model = RandomForestClassifier(random_state=42, n_jobs=-1)
|
||||||
|
|
||||||
|
# Tạo pipeline
|
||||||
|
pipeline = Pipeline([
|
||||||
|
# ('imputer', SimpleImputer(strategy='mean')),
|
||||||
|
('scaler', StandardScaler()),
|
||||||
|
('classifier', base_model),
|
||||||
|
])
|
||||||
|
# Thiết lập các tham số bạn muốn tối ưu hóa
|
||||||
|
param_grid = {
|
||||||
|
'classifier__n_estimators': [100, 300, 500, 700, 1000],
|
||||||
|
'classifier__max_depth': [6, 8, 10, 15, 20],
|
||||||
|
'classifier__criterion': ['gini', 'entropy'],
|
||||||
|
}
|
||||||
|
|
||||||
|
# Sử dụng GridSearchCV để tìm bộ tham số tốt nhất
|
||||||
|
grid_search = GridSearchCV(pipeline, param_grid, cv=5, scoring='accuracy', n_jobs=-1)
|
||||||
|
grid_search.fit(X_train, y_train)
|
||||||
|
|
||||||
|
# In ra bộ tham số tốt nhất
|
||||||
|
best_params = grid_search.best_params_
|
||||||
|
print("Best Parameters:", best_params)
|
||||||
|
|
||||||
|
# Dự đoán trên tập kiểm tra
|
||||||
|
y_pred = grid_search.predict(X_val)
|
||||||
|
|
||||||
|
# Đánh giá kết quả
|
||||||
|
accuracy = accuracy_score(y_val, y_pred)
|
||||||
|
print(f"Accuracy: {round(accuracy, 2)*100} %")
|
||||||
|
return grid_search
|
||||||
|
|
||||||
|
|
||||||
|
def save_model(name_file, grid_search):
|
||||||
|
dir_save_model = "model_train"
|
||||||
|
if not os.path.exists(dir_save_model):
|
||||||
|
os.mkdir(dir_save_model)
|
||||||
|
joblib.dump(grid_search, os.path.join(dir_save_model, name_file))
|
||||||
|
print("Done!")
|
||||||
|
|
||||||
|
|
||||||
|
def predict(model, data_crs, ndvi, vh, vv):
|
||||||
|
data_predict = []
|
||||||
|
for i in range(ndvi.shape[1]):
|
||||||
|
ndvi_tmp = ndvi.isel(y=i).values
|
||||||
|
vh_data = vh.sel(y=ndvi.y.values[i], method='nearest').values
|
||||||
|
vv_data = vv.sel(y=ndvi.y.values[i], method='nearest').values
|
||||||
|
all_tmp = np.concatenate((ndvi_tmp, vh_data, vv_data), axis=0)
|
||||||
|
data_predict.extend(all_tmp.T)
|
||||||
|
y_pred = model.predict(data_predict)
|
||||||
|
final_label = y_pred.reshape(ndvi.y.shape[0], ndvi.x.shape[0])
|
||||||
|
|
||||||
|
final_xarray_save = xr.DataArray(final_label, dims=("y", "x"))
|
||||||
|
final_xarray_save = final_xarray_save.rio.write_crs(data_crs)
|
||||||
|
|
||||||
|
x_values = ndvi.x.values
|
||||||
|
y_values = ndvi.y.values
|
||||||
|
|
||||||
|
data_array = xr.DataArray(final_xarray_save,
|
||||||
|
coords={'x': x_values, 'y': y_values},
|
||||||
|
dims=['y', 'x'])
|
||||||
|
data_array = data_array.rio.write_crs(ndvi.rio.crs)
|
||||||
|
return data_array
|
||||||
|
|
||||||
|
|
||||||
|
def cut_according_shp(thuanhoa_path, average_ndvi, data_array):
|
||||||
|
gdf = gpd.read_file(thuanhoa_path)
|
||||||
|
gdf = gdf.to_crs(average_ndvi.rio.crs)
|
||||||
|
polygon_coords = list(gdf.geometry.values[0].exterior.coords)
|
||||||
|
polygon_coordinates = [(x, y) for x, y in polygon_coords]
|
||||||
|
|
||||||
|
geometries = [
|
||||||
|
{
|
||||||
|
'type': 'Polygon',
|
||||||
|
'coordinates': [polygon_coordinates]
|
||||||
|
}
|
||||||
|
]
|
||||||
|
region_result = data_array.rio.clip(geometries, data_array.rio.crs, drop=False)
|
||||||
|
region_result = region_result.where(region_result >= 0, float('nan'))
|
||||||
|
return region_result
|
||||||
|
|
||||||
|
|
||||||
|
def compare(KD_path, KetQuaPhanLoaiDat, CODE_MAP, HT_MAP):
|
||||||
|
gdf = gpd.read_file(KD_path, crs="EPSG:9209")
|
||||||
|
polygon = gdf.geometry.values
|
||||||
|
label = gdf.tenchu.values
|
||||||
|
ouput_image = rioxarray.open_rasterio(KetQuaPhanLoaiDat)
|
||||||
|
code_tq = HT_MAP["TQ"]["data"][0]
|
||||||
|
code_pnn = HT_MAP["PNN"]["data"][0]
|
||||||
|
result = {}
|
||||||
|
for key, values in HT_MAP.items():
|
||||||
|
print(f"process {key}")
|
||||||
|
array_list = []
|
||||||
|
for i in range(len(polygon)):
|
||||||
|
po = polygon[i]
|
||||||
|
lb = label[i]
|
||||||
|
code_lb = CODE_MAP.get(lb, code_tq)
|
||||||
|
try:
|
||||||
|
qr = ouput_image.rio.clip([po], "EPSG:9209")
|
||||||
|
if code_lb in values["data"]:
|
||||||
|
if code_lb == code_pnn:
|
||||||
|
qr = qr.where((qr != float(code_pnn)), np.nan)
|
||||||
|
# qr = qr.where((qr != 3.0), np.nan)
|
||||||
|
elif code_lb == code_tq:
|
||||||
|
qr = qr.where((qr != float(code_pnn)), np.nan)
|
||||||
|
qr = qr.where((qr != 3.0), np.nan)
|
||||||
|
else:
|
||||||
|
qr = qr.where(qr != float(code_lb), np.nan)
|
||||||
|
else:
|
||||||
|
qr.values[:, :, :] = np.nan
|
||||||
|
array_list.append(qr)
|
||||||
|
except Exception as e:
|
||||||
|
pass
|
||||||
|
result.update({key: array_list})
|
||||||
|
return result
|
||||||
|
|
||||||
|
|
||||||
|
def save_result(result, HT_MAP):
|
||||||
|
# cmap = ListedColormap(colors)
|
||||||
|
save_path = "ThuanHoa/KetQua"
|
||||||
|
if not os.path.exists(save_path):
|
||||||
|
os.mkdir(save_path)
|
||||||
|
|
||||||
|
for k, v in result.items():
|
||||||
|
rs = merge_arrays(v, nodata = np.nan)
|
||||||
|
rs.rio.to_raster(f"{save_path}/{k}.tif")
|
||||||
|
print(f"save {save_path}/{k}.tif")
|
||||||
|
# img = rs.plot(cmap=cmap, add_colorbar=False)
|
||||||
|
# cbar = plt.colorbar(img)
|
||||||
|
# cbar.ax.set_yticklabels(labels)
|
||||||
|
# plt.title(f'{HT_MAP[k]["name"]}')
|
||||||
|
# plt.axis('off')
|
||||||
|
# plt.show()
|
||||||
|
|
||||||
|
|
||||||
|
def accuracy_test(test, data_array):
|
||||||
|
# cấu hình nhãn dữ liệu
|
||||||
|
label_mapping = {
|
||||||
|
"Lua tom": "0",
|
||||||
|
"Lua": "1",
|
||||||
|
"CHN": "2",
|
||||||
|
"CLN": "3",
|
||||||
|
"TS": "4",
|
||||||
|
"Song": "5",
|
||||||
|
"Dat xay dung": "6",
|
||||||
|
"Rung": "7"
|
||||||
|
}
|
||||||
|
|
||||||
|
chk = []
|
||||||
|
pred = []
|
||||||
|
dd = []
|
||||||
|
for idx, point in test.iterrows():
|
||||||
|
label = point.LULC
|
||||||
|
predict = data_array.sel(x=point.geometry.x, y=point.geometry.y, method='nearest').values
|
||||||
|
pred.append(label_mapping[label])
|
||||||
|
dd.append(str(predict))
|
||||||
|
chk.append(predict == int(label_mapping[label]))
|
||||||
|
test["code"] = pred
|
||||||
|
test["dd"] = dd
|
||||||
|
test["check"] = chk
|
||||||
|
path = "ThuanHoa/TestAccuracy"
|
||||||
|
if not os.path.exists(path):
|
||||||
|
os.mkdir(path)
|
||||||
|
test.to_file(f"{path}/result.shp")
|
||||||
|
|
||||||
|
percentage_true = np.mean(chk) * 100
|
||||||
|
print(f"độ chính xác: {percentage_true:.2f}%")
|
||||||
+1
Submodule backup_code_training/CSIROBoeingPhase5-Vietnam added at 8f0cb55cba
@@ -9,7 +9,6 @@ from typing import Tuple, Optional, Dict
|
|||||||
from sklearn.neighbors import KNeighborsRegressor
|
from sklearn.neighbors import KNeighborsRegressor
|
||||||
from sklearn.ensemble import RandomForestRegressor
|
from sklearn.ensemble import RandomForestRegressor
|
||||||
import warnings
|
import warnings
|
||||||
from pathlib import Path
|
|
||||||
warnings.filterwarnings('ignore')
|
warnings.filterwarnings('ignore')
|
||||||
|
|
||||||
|
|
||||||
@@ -376,22 +375,6 @@ class DeepInpaintingStrategy(CloudRemovalStrategy):
|
|||||||
# Convert to tensor and add batch dimension
|
# Convert to tensor and add batch dimension
|
||||||
input_tensor = torch.from_numpy(input_array).unsqueeze(0).to(self.device)
|
input_tensor = torch.from_numpy(input_array).unsqueeze(0).to(self.device)
|
||||||
|
|
||||||
if hasattr(self.model, 'encoder'):
|
|
||||||
# Custom UNet from train_cloud_removal.py
|
|
||||||
expected_channels = self.model.encoder[0].double_conv[0].in_channels
|
|
||||||
elif hasattr(self.model, 'inc') and hasattr(self.model.inc.double_conv[0], 'in_channels'):
|
|
||||||
expected_channels = self.model.inc.double_conv[0].in_channels
|
|
||||||
elif hasattr(self.model, 'conv1') and hasattr(self.model.conv1, 'in_channels'):
|
|
||||||
expected_channels = self.model.conv1.in_channels
|
|
||||||
else:
|
|
||||||
expected_channels = 6
|
|
||||||
|
|
||||||
|
|
||||||
if expected_channels > input_tensor.shape[1]:
|
|
||||||
pad_channels = expected_channels - input_tensor.shape[1]
|
|
||||||
padding = torch.zeros(1, pad_channels, *input_tensor.shape[2:]).to(self.device)
|
|
||||||
input_tensor = torch.cat([input_tensor, padding], dim=1)
|
|
||||||
|
|
||||||
# Run through U-Net
|
# Run through U-Net
|
||||||
with torch.no_grad():
|
with torch.no_grad():
|
||||||
output_tensor = self.model(input_tensor)
|
output_tensor = self.model(input_tensor)
|
||||||
@@ -1,93 +0,0 @@
|
|||||||
import os
|
|
||||||
import torch
|
|
||||||
import numpy as np
|
|
||||||
import xarray as xr
|
|
||||||
from torch.utils.data import Dataset
|
|
||||||
import glob
|
|
||||||
|
|
||||||
class NDVITimeSeriesDataset(Dataset):
|
|
||||||
def __init__(self, sequence_length=3, spatial=False):
|
|
||||||
"""
|
|
||||||
Đọc dữ liệu S2 từ cache, tính NDVI và tạo Time-Series.
|
|
||||||
spatial=False -> Output 1D cho LSTM/ARIMA
|
|
||||||
spatial=True -> Output 2D cho ConvLSTM
|
|
||||||
"""
|
|
||||||
self.sequence_length = sequence_length
|
|
||||||
self.spatial = spatial
|
|
||||||
self.data_seqs = []
|
|
||||||
self.targets = []
|
|
||||||
|
|
||||||
# Load from cache
|
|
||||||
cache_files = glob.glob("dataset_cache/*.nc")
|
|
||||||
s2_files = [f for f in cache_files if len(os.path.basename(f)) == 35] # S2 cache filenames usually have length 32 + 3 (.nc)
|
|
||||||
|
|
||||||
if not s2_files:
|
|
||||||
print("[WARNING] Không tìm thấy dữ liệu S2 trong cache! Dùng dummy data.")
|
|
||||||
self._create_dummy()
|
|
||||||
return
|
|
||||||
|
|
||||||
try:
|
|
||||||
print(f"[DATA] Loading real data from {s2_files[0]}")
|
|
||||||
ds = xr.open_dataset(s2_files[0], engine='netcdf4')
|
|
||||||
if 'time' not in ds.dims or len(ds.time) < sequence_length + 1:
|
|
||||||
self._create_dummy()
|
|
||||||
return
|
|
||||||
|
|
||||||
# Tính NDVI: (B08 - B04) / (B08 + B04)
|
|
||||||
b8 = ds['B08'].astype(np.float32)
|
|
||||||
b4 = ds['B04'].astype(np.float32)
|
|
||||||
ndvi = (b8 - b4) / (b8 + b4 + 1e-8)
|
|
||||||
ndvi = ndvi.fillna(0).values # shape: (time, y, x)
|
|
||||||
|
|
||||||
# Lấy 1 pixel trung tâm hoặc toàn bộ ảnh
|
|
||||||
if not self.spatial:
|
|
||||||
# Average pooling over space for 1D time series
|
|
||||||
ndvi = ndvi.mean(axis=(1, 2)) # shape: (time,)
|
|
||||||
for i in range(len(ndvi) - sequence_length):
|
|
||||||
self.data_seqs.append(ndvi[i:i+sequence_length])
|
|
||||||
self.targets.append(ndvi[i+sequence_length])
|
|
||||||
else:
|
|
||||||
# Spatial data for ConvLSTM
|
|
||||||
# Downsample to 64x64 to avoid OOM
|
|
||||||
from skimage.transform import resize
|
|
||||||
T = len(ndvi)
|
|
||||||
ndvi_resized = np.zeros((T, 64, 64))
|
|
||||||
for t in range(T):
|
|
||||||
ndvi_resized[t] = resize(ndvi[t], (64, 64))
|
|
||||||
|
|
||||||
for i in range(T - sequence_length):
|
|
||||||
self.data_seqs.append(ndvi_resized[i:i+sequence_length]) # (seq, 64, 64)
|
|
||||||
self.targets.append(ndvi_resized[i+sequence_length]) # (64, 64)
|
|
||||||
|
|
||||||
except Exception as e:
|
|
||||||
print(f"[ERROR] {e}. Dùng dummy data.")
|
|
||||||
self._create_dummy()
|
|
||||||
|
|
||||||
def _create_dummy(self):
|
|
||||||
T = 20
|
|
||||||
if not self.spatial:
|
|
||||||
ndvi = np.random.rand(T).astype(np.float32)
|
|
||||||
for i in range(T - self.sequence_length):
|
|
||||||
self.data_seqs.append(ndvi[i:i+self.sequence_length])
|
|
||||||
self.targets.append(ndvi[i+self.sequence_length])
|
|
||||||
else:
|
|
||||||
ndvi = np.random.rand(T, 64, 64).astype(np.float32)
|
|
||||||
for i in range(T - self.sequence_length):
|
|
||||||
self.data_seqs.append(ndvi[i:i+self.sequence_length])
|
|
||||||
self.targets.append(ndvi[i+self.sequence_length])
|
|
||||||
|
|
||||||
def __len__(self):
|
|
||||||
return len(self.data_seqs)
|
|
||||||
|
|
||||||
def __getitem__(self, idx):
|
|
||||||
x = torch.tensor(self.data_seqs[idx], dtype=torch.float32)
|
|
||||||
y = torch.tensor(self.targets[idx], dtype=torch.float32)
|
|
||||||
|
|
||||||
if not self.spatial:
|
|
||||||
x = x.unsqueeze(1) # (seq_len, features=1)
|
|
||||||
y = y.unsqueeze(0) # (1,)
|
|
||||||
else:
|
|
||||||
x = x.unsqueeze(1) # (seq_len, channels=1, H, W)
|
|
||||||
y = y.unsqueeze(0) # (1, H, W)
|
|
||||||
|
|
||||||
return x, y
|
|
||||||
Binary file not shown.
Binary file not shown.
Some files were not shown because too many files have changed in this diff Show More
Reference in New Issue
Block a user