البرمجة الجيومكانية

import arcpy
arcpy.env.overwriteOutput = True

# نصيحة عملية:
# ثبّت workspace واضح، واستخدم geodatabases بدل مسارات مبعثرة
arcpy.env.workspace = r"C:\GIS\Project\work.gdb"

project/
  data_raw/        # بيانات خام (لا تعدّلها)
  data_processed/  # بيانات نظيفة/موحّدة
  scripts/         # سكربتات
  notebooks/       # تجارب وتحليل
  outputs/         # خرائط، جداول، تقارير
  docs/            # توثيق ومنهجية
  README.md        # وصف المشروع وكيفية التشغيل

import arcpy
arcpy.env.overwriteOutput = True

study = r"C:\GIS\work.gdb\boundary"
roads = r"C:\GIS\data.gdb\roads"
out_gdb = r"C:\GIS\out.gdb"

clip_out = arcpy.analysis.Clip(roads, study, f"{out_gdb}\\roads_clip")
buf_out  = arcpy.analysis.Buffer(clip_out, f"{out_gdb}\\roads_buf_50m", "50 Meters")

arcpy.management.AddField(clip_out, "len_m", "DOUBLE")
arcpy.management.CalculateGeometryAttributes(
    clip_out, [["len_m","LENGTH_GEODESIC"]], length_unit="METERS"
)

print("Completed.")

import rasterio
import numpy as np

with rasterio.open("data_processed/dem.tif") as src:
    arr = src.read(1)
    print("CRS:", src.crs)
    print("Transform:", src.transform)
    print("Min/Max:", float(np.nanmin(arr)), float(np.nanmax(arr)))

import rasterio
import numpy as np

nir_path = "data_processed/S2_B08.tif"  # NIR
red_path = "data_processed/S2_B04.tif"  # Red

with rasterio.open(nir_path) as nir, rasterio.open(red_path) as red:
    N = nir.read(1).astype("float32")
    R = red.read(1).astype("float32")

    ndvi = (N - R) / (N + R + 1e-6)

    profile = nir.profile
    profile.update(dtype="float32", count=1)

    with rasterio.open("outputs/ndvi.tif", "w", **profile) as dst:
        dst.write(ndvi, 1)

import rioxarray as rxr
import numpy as np

dem = rxr.open_rasterio("data_processed/dem.tif").squeeze()

dx = dem.differentiate("x")
dy = dem.differentiate("y")

slope = np.degrees(np.arctan(np.hypot(dx, dy)))
out = dem.copy(data=slope)
out.rio.to_raster("outputs/slope_deg.tif")

import geopandas as gpd
import pandas as pd

districts = gpd.read_file("data_processed/districts.gpkg").to_crs(epsg=32638)
roads = gpd.read_file("data_processed/roads.gpkg", layer="roads").to_crs(districts.crs)

# انضمام مكاني: أي طريق يقع ضمن أي حي
join = gpd.sjoin(roads[["geometry"]], districts[["name","geometry"]], how="left", predicate="intersects")
summary = (join.dropna(subset=["name"])
              .assign(len_m = join.geometry.length)
              .groupby("name")["len_m"].sum()
              .sort_values(ascending=False)
              .reset_index())

summary["len_km"] = (summary["len_m"]/1000).round(3)
summary.drop(columns=["len_m"]).to_csv("outputs/roads_length_by_district.csv", index=False)

<link rel="stylesheet" href="https://unpkg.com/leaflet@1.9.4/dist/leaflet.css">
<div id="map" style="height:420px;border-radius:14px"></div>
<script src="https://unpkg.com/leaflet@1.9.4/dist/leaflet.js"></script>
<script>
  const map = L.map('map').setView([24.7,46.7], 6);

  const osm = L.tileLayer('https://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png', {
    maxZoom: 19, attribution:'© OpenStreetMap'
  }).addTo(map);

  fetch('data/points.geojson')
    .then(r => r.json())
    .then(g => {
      L.geoJSON(g, {
        onEachFeature: (f, layer) => {
          const name = (f.properties && f.properties.name) ? f.properties.name : '';
          layer.bindPopup(name);
        }
      }).addTo(map);
    });
</script>

<link rel="stylesheet" href="https://js.arcgis.com/4.29/esri/themes/light/main.css">
<div id="view" style="height:420px"></div>
<script src="https://js.arcgis.com/4.29/"></script>
<script>
require([
  "esri/Map",
  "esri/views/MapView",
  "esri/layers/FeatureLayer"
], (Map, MapView, FeatureLayer) => {

  const map = new Map({ basemap: "streets-navigation-vector" });

  const view = new MapView({
    container: "view",
    map,
    center: [46.7, 24.7],
    zoom: 6
  });

  const layer = new FeatureLayer({
    url: "YOUR_FEATURESERVER_URL/0"
  });

  map.add(layer);
});
</script>

CREATE DATABASE gisdb;
\c gisdb
CREATE EXTENSION postgis;

SELECT postgis_full_version();

CREATE INDEX idx_roads_geom
ON public.roads
USING GIST (geom);

SELECT d.id, d.name,
       SUM(ST_Length(ST_Intersection(r.geom, d.geom))) / 1000 AS km_len
FROM districts d
JOIN roads r
  ON ST_Intersects(r.geom, d.geom)
GROUP BY d.id, d.name
ORDER BY km_len DESC;

// مثال WMS
const wms = L.tileLayer.wms("https://YOUR-GEOSERVER/geoserver/wms", {
  layers: "workspace:layername",
  format: "image/png",
  transparent: true
}).addTo(map);