New algo

backend/src/main.rs

@@ -17,7 +17,6 @@ async fn main() {
     // Load DEM dataset and generate the base image
     let ms = MapSegment::load_from_dataset("USGS_13_n44w072_20240617.tif")
         .expect("Failed to load dataset");
-    //ms.gen_map().expect("Failed to generate map");
 
     let ms = Arc::new(RwLock::new(ms));
 

backend/src/map_segment.rs

@@ -13,6 +13,7 @@ pub struct MapSegment {
     pub height: usize,
     pub height_map: Buffer<f32>,
     pub img: RgbaImage,
+    pub los_map: Option<Vec<bool>>,
 }
 
 impl MapSegment {
@@ -44,7 +45,8 @@ impl MapSegment {
             width,
             height,
             height_map,
-            img
+            img,
+            los_map: None,
         })
     }
 
@@ -242,7 +244,6 @@ impl MapSegment {
         image::imageops::resize(&cropped, tile_size, tile_size, FilterType::Lanczos3)
     }
 
-#####
     pub fn has_line_of_sight(
         &self,
         from_lat: f64, // Latitude of the observer
@@ -277,7 +278,53 @@ impl MapSegment {
         true // Line of sight is clear
     }
 
-    fn latlon_to_pixel(&self, lat: f64, lon: f64) -> (usize, usize) {
+    /// Efficiently generate a line-of-sight map from the observer (in lat/lon) out to a given radius (in pixels)
+    /// Now accounts for Earth's curvature.
+    pub fn generate_los_map(&mut self, observer_lat: f64, observer_lon: f64, radius: usize, observer_height: f32) {
+        let mut los_map = vec![false; self.width * self.height];
+        let (obs_px, obs_py) = self.latlon_to_pixel(observer_lat, observer_lon);
+        let obs_elev = self.get_elevation(obs_px, obs_py).unwrap_or(0.0) + observer_height;
+
+        // Constants for Earth's curvature
+        const EARTH_RADIUS_M: f64 = 6_371_000.0; // meters
+        let pixel_size_m = self.pixel_width.hypot(self.pixel_height).abs(); // crude average pixel size in meters
+
+        // Cast rays in all directions from the observer
+        let num_rays = 360.max(radius * 8); // More rays for larger radius
+        for angle_step in 0..num_rays {
+            let theta = (angle_step as f64) * std::f64::consts::TAU / (num_rays as f64);
+            let dx = theta.cos();
+            let dy = theta.sin();
+            let mut max_angle = std::f32::NEG_INFINITY;
+            for r in 1..=radius {
+                let xi = obs_px as f64 + dx * r as f64;
+                let yi = obs_py as f64 + dy * r as f64;
+                let x = xi.round() as isize;
+                let y = yi.round() as isize;
+                if x < 0 || y < 0 || x >= self.width as isize || y >= self.height as isize {
+                    break;
+                }
+                let x = x as usize;
+                let y = y as usize;
+                let idx = y * self.width + x;
+                let elev = self.get_elevation(x, y).unwrap_or(0.0);
+                let dist_m = (r as f64) * pixel_size_m;
+                // Drop due to Earth's curvature (in meters)
+                let curvature_drop = dist_m * dist_m / (2.0 * EARTH_RADIUS_M);
+                // Adjusted elevation at this point
+                let adj_elev = elev - curvature_drop as f32;
+                let dist = dist_m as f32;
+                let angle = if dist > 0.0 { (adj_elev - obs_elev) / dist } else { std::f32::NEG_INFINITY };
+                if angle > max_angle {
+                    los_map[idx] = true;
+                    max_angle = angle;
+                }
+            }
+        }
+        self.los_map = Some(los_map);
+    }
+
+    pub fn latlon_to_pixel(&self, lat: f64, lon: f64) -> (usize, usize) {
         let px = ((lon - self.origin_x) / self.pixel_width) as usize;
         let py = ((lat - self.origin_y) / self.pixel_height) as usize;
 
@@ -314,7 +361,7 @@ fn tile_to_latlon(x: u32, y: u32, zoom: u8) -> (f64, f64) {
     (lat_deg, lon_deg)
 }
 
-fn tile_bbox_latlon(x: u32, y: u32, zoom: u8) -> ((f64, f64), (f64, f64)) {
+pub fn tile_bbox_latlon(x: u32, y: u32, zoom: u8) -> ((f64, f64), (f64, f64)) {
     // Top-left corner
     let (lat1, lon1) = tile_to_latlon(x, y, zoom);
     // Bottom-right corner (x+1, y+1)

backend/src/web/get_tile_los.rs

@@ -1,7 +1,6 @@
 use axum::{extract::{Path, State}, http::StatusCode, response::{IntoResponse, Response}};
 use std::{sync::{Arc, RwLock}, io::Cursor};
 use crate::map_segment::MapSegment;
-use gdal::raster::Buffer;
 use image::{ImageFormat, RgbaImage, ImageBuffer, Rgba, DynamicImage};
 use bytes::Bytes;
 
@@ -43,9 +42,4 @@ pub async fn handler(
         ],
         bytes
     ).into_response()
-}
-
-fn within_banded_range(n: f32, multiple: i32, buffer: f64) -> bool {
-    let remainder = ( n as f64 ) % multiple as f64;
-    remainder <= buffer || remainder >= ( multiple as f64 - buffer )
-}
+}

backend/src/web/get_tile_v2.rs

@@ -1,7 +1,6 @@
 use axum::{extract::{Path, State}, http::StatusCode, response::{IntoResponse, Response}};
 use std::{sync::{Arc, RwLock}, io::Cursor};
 use crate::map_segment::MapSegment;
-use gdal::raster::Buffer;
 use image::{ImageFormat, RgbaImage, ImageBuffer, Rgba, DynamicImage};
 use bytes::Bytes;
 

backend/src/web/get_tile_v3.rs

@@ -0,0 +1,62 @@
+use axum::{extract::{Path, State}, http::StatusCode, response::{IntoResponse, Response}};
+use std::{sync::{Arc, RwLock}, io::Cursor};
+use crate::map_segment::MapSegment;
+use image::{ImageFormat, RgbaImage, ImageBuffer, Rgba, DynamicImage};
+use bytes::Bytes;
+
+pub async fn handler(
+    Path((z, x, y)): Path<(u8, u32, u32)>,
+    State(ms): State<Arc<RwLock<MapSegment>>>,
+) -> Response {
+    let ms = match ms.read() {
+        Ok(m) => m,
+        Err(_) => return StatusCode::INTERNAL_SERVER_ERROR.into_response(),
+    };
+
+    let mut img: RgbaImage = ImageBuffer::from_pixel(256, 256, Rgba([255, 255, 255, 0]));
+
+    // Only proceed if los_map is present
+    if let Some(los_map) = &ms.los_map {
+        // Calculate tile bounds in raster pixel space
+        let ((lat_north, lon_west), (lat_south, lon_east)) = crate::map_segment::tile_bbox_latlon(x, y, z);
+        let (px_min, py_min) = ms.latlon_to_pixel(lat_north, lon_west);
+        let (px_max, py_max) = ms.latlon_to_pixel(lat_south, lon_east);
+        let step_size_x: f32 = (px_max - px_min) as f32 / 256.0;
+        let step_size_y: f32 = (py_max - py_min) as f32 / 256.0;
+        let px_min_f32 = px_min as f32;
+        let py_min_f32 = py_min as f32;
+
+        for x in 0..256 {
+            for y in 0..256 {
+                let source_x = px_min_f32 + (x as f32 * step_size_x);
+                let source_y = py_min_f32 + (y as f32 * step_size_y);
+                let px = source_x.round() as usize;
+                let py = source_y.round() as usize;
+                if px < ms.width && py < ms.height {
+                    let idx = py * ms.width + px;
+                    if idx < los_map.len() && los_map[idx] {
+                        img.put_pixel(x as u32, y as u32, Rgba([190, 0, 150, 255]));
+                    }
+                }
+            }
+        }
+    }
+
+    let mut buf = Cursor::new(Vec::new());
+    if DynamicImage::ImageRgba8(img)
+        .write_to(&mut buf, ImageFormat::Png)
+        .is_err()
+    {
+        println!("Failed to encode tile PNG");
+        return StatusCode::INTERNAL_SERVER_ERROR.into_response();
+    }
+
+    let bytes = Bytes::from(buf.into_inner());
+    (
+        [
+            ("Content-Type", "image/png"),
+            ("Content-Length", &bytes.len().to_string()),
+        ],
+        bytes
+    ).into_response()
+}

backend/src/web/mod.rs

@@ -3,6 +3,7 @@ pub mod run_scan;
 pub mod get_tile;
 pub mod get_tile_v2;
 pub mod get_tile_los;
+pub mod get_tile_v3;
 
 use std::sync::{Arc, RwLock};
 use crate::map_segment::MapSegment;
@@ -16,6 +17,7 @@ pub fn router(ms: Arc<RwLock<MapSegment>>) -> Router {
         .route("/tiles_rf/{z}/{x}/{y}/tile.png", get(get_tile::handler))
         .route("/tiles_los/{z}/{x}/{y}/tile.png", get(get_tile_los::handler))
         .route("/tiles/{z}/{x}/{y}/tile.png", get(get_tile_v2::handler))
+        .route("/tile_v3/{z}/{x}/{y}/tile.png", get(get_tile_v3::handler))
         .route("/los", post(get_los::handler))
         .route("/scan", post(run_scan::handler))
         .route("/", axum::routing::get(|| async { "Hello, World!" }))

backend/src/web/run_scan.rs

@@ -9,37 +9,20 @@ use crate::map_segment::MapSegment;
 pub struct ScanRequest {
     from_lat: f64,
     from_lon: f64,
-    polygon: Vec<(f64, f64)>,
+    radius: usize,
+    observer_height: f32,
 }
 
 pub async fn handler(
     State(ms): State<Arc<RwLock<MapSegment>>>,
     Json(payload): Json<ScanRequest>,
 ) -> Result<Json<bool>, axum::http::StatusCode> {
-    println!("Starting scan...");
+    println!("Starting LOS map generation...");
     let mut ms = ms.write().map_err(|_| axum::http::StatusCode::INTERNAL_SERVER_ERROR)?;
 
-    let pixels = ms.scan_polygon(&payload.polygon);
-    let mut count = 0;
+    ms.generate_los_map(payload.from_lat, payload.from_lon, payload.radius, payload.observer_height);
 
-    for (x, y, _elev) in pixels {
-        let (to_lat, to_lon) = ms.pixel_to_latlon(x, y);
-        let has_los = ms.has_line_of_sight(
-            payload.from_lat,
-            payload.from_lon,
-            to_lat,
-            to_lon,
-            5.0,
-            10.0,
-        );
-
-        if has_los {
-            ms.draw_pixel(x as u32, y as u32);
-            count += 1;
-        }
-    }
-
-    println!("Scan complete: {} pixels with LOS", count);
+    println!("LOS map generated for observer at ({}, {}) radius {} height {}", payload.from_lat, payload.from_lon, payload.radius, payload.observer_height);
     Ok(Json(true))
 }
 

frontend/src/App.vue

@@ -9,12 +9,20 @@
 
       <l-marker v-for="marker, index in markers" :lat-lng="marker" @click="removeMarker(index)"></l-marker>
 
+      <l-tile-layer
+        url="http://localhost:3000/tile_v3/{z}/{x}/{y}/tile.png"
+        layer-type="base"
+        opacity="0.45"
+        name="Elevation"
+      ></l-tile-layer>
+<!--
       <l-tile-layer
         url="http://localhost:3000/tiles_los/{z}/{x}/{y}/tile.png"
         layer-type="base"
         opacity="0.45"
         name="Elevation"
       ></l-tile-layer>
+-->
     </l-map>
   </div>
   <button @click="selection" >Select Area</button>
@@ -87,7 +95,8 @@ export default {
       let res = await axios.post('http://localhost:3000/scan', {
         from_lat: this.markers[0].lat,
         from_lon: this.markers[0].lng,
-        polygon: this.polyPoints.map((p) => [p.lat, p.lng])
+        radius: 4000,
+        observer_height: 10
       });
 
       console.log(res);