diff --git a/README.md b/README.md
index 18fb5b5..1010ca1 100644
--- a/README.md
+++ b/README.md
@@ -6,10 +6,10 @@
 [![Build](https://github.com/prina404/raycast2D/actions/workflows/build.yml/badge.svg)](https://github.com/prina404/raycast2D/actions/workflows/build.yml)
 [![Tests](https://github.com/prina404/raycast2D/actions/workflows/tests.yml/badge.svg)](https://github.com/prina404/raycast2D/actions/workflows/tests.yml)
 
-`raycast2D` is a fast, single-core 2D raycasting implementation written in C with a small Python API.
+`raycast2D` is a fast, single-thread 2D raycasting implementation written in C with a small Python API.
 It operates on NumPy occupancy grids / binary images (free cells are non-zero; occupied cells are `0`) and computes ray intersections using Bresenham line algorithm.
 
-The core package depends only on `numpy`. Optional extras are provided for the interactive demo and development.
+The core package depends only on `numpy`.
 
 ![](media/raycast_compressed.gif)
 
@@ -18,30 +18,66 @@ The core package depends only on `numpy`. Optional extras are provided for the i
 ```bash
 $ pip install raycast2D
 ```
-### Basic raycast on an image
+### Basic `cast` function usage
 
 ```python
 import numpy as np
-from raycast2D import cast
 from PIL import Image
-import matplotlib.pyplot as plt
+from raycast2D import cast
 
 img = Image.open("<path/to/img.png>")
-img_array = np.array(img)
-img_array = img_array[:, :, 0].astype(np.uint8)  # Use single channel
-# img_array = img_array[img_array < 100] = 0  # threshold obstacles if needed
+img_array = np.array(img)[:, :, 0] # Use single channel
+pose = (250, 250, 0.5)  # x, y, yaw in radians 
+
+rays = cast(img_array, pose, num_rays=360, ray_length=500)
+
+print(rays[:5])
+# Expected output (example):
+# [[x0 y0]
+#  [x1 y1]
+#  [x2 y2]
+#  [x3 y3]
+#  [x4 y4]]
+```
+
+### Using the `Lidar2D` class
 
-rays = cast(img_array, pose=(250, 250), num_rays=360, ray_length=500)
-# rays is an array of shape (N, 2) where each row contains the (x, y) coordinates of the ray collisions
+```python
+from raycast2D import Lidar2D
+# [...] 
+lidar = Lidar2D(num_rays=360, FOV=360, ray_length=500)
+lidar.set_map(img_array)
 
-plt.imshow(img_array, cmap='gray')
-plt.scatter([250], [250], c='green', s=10)
-plt.scatter(rays[:, 0], rays[:, 1], c='blue', s=1)
+rays = lidar.scan((250, 250, 0.5))
+
+print(rays[:5])
+# same output as before...
+```
+
+### Plotting example
+
+```python
+import numpy as np
+import matplotlib.pyplot as plt
+from PIL import Image
+from raycast2D import cast
+
+img = Image.open("media/lab_intel.png")
+img_array = np.array(img)[:, :, 0] # Use single channel
+
+pose = (350, 350, -1.0)
+rays = cast(img_array, pose=pose, num_rays=200, FOV=90, ray_length=250)
+
+plt.imshow(img_array, cmap="gray")
+plt.scatter([pose[0]], [pose[1]], c="green", s=10)
+plt.scatter(rays[:, 0], rays[:, 1], c="blue", s=1)
 for ray in rays:
-    plt.plot([250, ray[0]], [250, ray[1]], c='red', linewidth=0.5, alpha=0.3)
+    plt.plot([pose[0], ray[0]], [pose[1], ray[1]], c="red", linewidth=0.5, alpha=0.3)
 plt.show()
 ```
 
+![](media/example.png)
+
 ## Performance
 
 The benchmark script in `test/benchmark.py` runs a small set of examples and prints throughput in rays/s.
@@ -49,9 +85,9 @@ Example results (tested on an i7-9700k):
 
 |  Map size | Ray length | Rays per call | Mean time (ms) | Throughput (rays/s) |
 | --------: | ---------: | ------------: | -------------: | ------------------: |
-|   512×512 |       2000px |          5000 |         0.2904 |          17,216,880 |
-| 4096×4096 |       2000px |          5000 |         1.7626 |           2,836,661 |
-| 8192×8192 |       2000px |          5000 |         7.3452 |             680,719 |
+|   512×512 |     2000px |          5000 |         0.2490 |          20,083,924 |
+| 4096×4096 |     2000px |          5000 |         0.2848 |          17,556,549 |
+| 8192×8192 |     2000px |          5000 |         0.2571 |          19,446,324 |
 
 To reproduce on your machine:
 
@@ -59,8 +95,6 @@ To reproduce on your machine:
 $ python3 test/benchmark.py
 ```
 
-
-
 ## Interactive demo
 
 The repository includes an interactive `pygame` demo that raycasts from the current mouse position.
diff --git a/media/example.png b/media/example.png
new file mode 100644
index 0000000..fc65f03
Binary files /dev/null and b/media/example.png differ
diff --git a/src/raycast2D/__init__.py b/src/raycast2D/__init__.py
index 2e16ec8..70c2afd 100644
--- a/src/raycast2D/__init__.py
+++ b/src/raycast2D/__init__.py
@@ -1,3 +1,3 @@
-from .raycast2D import cast
+from .raycast2D import cast, Lidar2D
 
-__all__ = ["cast"]
\ No newline at end of file
+__all__ = ["cast", "Lidar2D"]
\ No newline at end of file
diff --git a/src/raycast2D/raycast2D.py b/src/raycast2D/raycast2D.py
index 2dd0024..7ddbc46 100644
--- a/src/raycast2D/raycast2D.py
+++ b/src/raycast2D/raycast2D.py
@@ -1,3 +1,4 @@
+from typing import Optional
 from . import raycaster
 import numpy as np
 from numpy.typing import NDArray
@@ -6,6 +7,94 @@
 import warnings
 
 
+class Lidar2D:
+
+    def __init__(
+        self,
+        num_rays: int = 1000,
+        FOV: int = 360,
+        ray_length: int = 500,
+        occupancy_grid: NDArray | None = None,
+        obstacle_threshold: int = 0,
+    ) -> None:
+        """
+        Initialize the Lidar2D raycaster.
+        Args:
+            num_rays: Number of rays to cast (samples uniformly over the FOV).
+            FOV: Field of view in degrees.
+            ray_length: Max ray length in pixels.
+            occupancy_grid: Optional occupancy grid to set at initialization.
+            obstacle_threshold: Cells with values less than or equal to this are treated as occupied.
+        """
+        self.angles = np.linspace(
+            -math.radians(FOV) / 2.0,
+            math.radians(FOV) / 2.0,
+            num_rays,
+        )
+        self.ray_length = ray_length
+        self.obstacle_threshold = obstacle_threshold
+        self._map = None
+        if occupancy_grid is not None:
+            self.set_map(occupancy_grid)
+
+    def scan(
+        self,
+        pose: tuple[int, int] | tuple[int, int, float],
+        image: Optional[NDArray] = None,
+        only_true_collisions: bool = True
+    ) -> NDArray[np.uint32]:
+        """
+        Perform a LiDAR scan from the given pose on the occupancy image.
+
+        Args:
+            pose: ``(x, y)`` or ``(x, y, yaw)`` where ``yaw`` is in radians.
+            image: Optional occupancy grid of shape ``(H, W)`` or ``(H, W, 1)`` with an integer dtype.
+                The raycaster treats values less than or equal to ``obstacle_threshold`` as occupied;
+                the start cell at ``pose`` must be free.
+                NOTE: the provided image overrides any occupancy grid set at initialization.
+            only_true_collisions: If True, return only rays that actually hit an obstacle
+                (LiDAR-style). If False, return endpoints for all rays.
+
+        Returns:
+            An array of shape ``(N, 2)`` with dtype ``uint32`` containing ``(x, y)`` collision coordinates.
+        """
+        if self._map is None and image is None:
+            raise ValueError(
+                "No occupancy grid provided at initialization or scan time.")
+            
+        if image is not None :
+            _map = (_check_map(image) > self.obstacle_threshold).astype(np.uint8)
+        else:
+            _map = self._map
+            
+        _check_pose(_map, pose)
+
+        p = _POSE(*pose)  # type: ignore
+        x_rays = np.round(p.x + self.ray_length *
+                          np.cos(self.angles + p.yaw)).astype(np.uint32)
+        y_rays = np.round(p.y + self.ray_length *
+                          np.sin(self.angles + p.yaw)).astype(np.uint32)
+
+        rays = np.column_stack((x_rays, y_rays, np.zeros_like(y_rays)))
+
+        # raycaster modifies the rays array inplace to include collisions
+        raycaster.raycast(_map, rays, p.x, p.y)
+
+        if only_true_collisions:    # LiDaR-style output
+            return rays[rays[:, -1] == 1][:, :-1]
+        return rays[:, :-1]
+
+    def set_map(self, occupancy_grid: NDArray) -> None:
+        """
+        Set or update the occupancy grid used by the raycaster.
+
+        Args:
+            occupancy_grid: Occupancy grid of shape ``(H, W)`` or ``(H, W, 1)`` with an integer dtype.
+        """
+        self._map = _check_map(occupancy_grid)
+        self._map = (self._map > self.obstacle_threshold).astype(np.uint8)
+
+
 def cast(
     image: NDArray,
     pose: tuple[int, int] | tuple[int, int, float],
@@ -31,14 +120,10 @@ def cast(
         An array of shape ``(N, 2)`` with dtype ``uint32`` containing ``(x, y)`` collision coordinates.
         If ``only_true_collisions=True``, ``N`` is the number of hits; otherwise ``N == num_rays``.
     """
-    _check_args(image, pose)
-    if len(image.shape) == 3:   # take only the first channel
-        image = image[:, :, 0]
-    if image.dtype != np.uint8:
-        warnings.warn(f"Converting image of dtype {image.dtype} to uint8 for raycasting.")
-        image = image.astype(np.uint8)
-
-    p = _POSE(*pose) # type: ignore
+    image = _check_map(image)
+    _check_pose(image, pose)
+
+    p = _POSE(*pose)  # type: ignore
     half_fov = math.radians(FOV) / 2.0
     start_angle = p.yaw - half_fov
     end_angle = p.yaw + half_fov
@@ -46,47 +131,61 @@ def cast(
     theta = np.linspace(start_angle, end_angle, num_rays)
     x_rays = np.round(p.x + ray_length * np.cos(theta)).astype(np.uint32)
     y_rays = np.round(p.y + ray_length * np.sin(theta)).astype(np.uint32)
-    
+
     rays = np.column_stack((x_rays, y_rays, np.zeros_like(y_rays)))
-    
+
     # raycaster modifies the rays array inplace to include collisions
     raycaster.raycast(image, rays, p.x, p.y)
-    
+
     if only_true_collisions:    # LiDaR-style output
-        return rays[rays[:, 2] == 1][:, :2]
-    return rays[:, :2]
+        return rays[rays[:, -1] == 1][:, :-1]
+    return rays[:, :-1]
 
 
 @dataclass
 class _POSE:
     x: int
     y: int
-    yaw: float = -math.pi / 2
-
-
-def _check_args(img: NDArray, pose: tuple) -> None:
-    if not np.issubdtype(img.dtype, np.integer):
-        raise ValueError(
-            f"Received an array with dtype {img.dtype}, The only supported dtypes are subtypes of np.integral.")
+    yaw: float = 0.0
+    
+    def __post_init__(self):
+        self.yaw -= math.pi / 2.0  # adjust so yaw=0 points upwards
 
-    shape = img.shape
-    if len(shape) < 2 or (len(shape) > 2 and shape[2] > 1):
-        raise ValueError(
-            f"Received an array of shape {shape}. The only supported shapes are (H x W), (H x W x 1)."
-        )
+def _check_pose(img: NDArray[np.uint8], pose: tuple) -> None:
     if len(pose) < 2 or len(pose) > 3:
         raise ValueError(
             f"Received a pose tuple {pose}. The only supported pose formats are (x, y) and (x, y, yaw).")
 
     x, y = pose[0], pose[1]
-    H, W = shape[0], shape[1]
+    H, W = img.shape[0], img.shape[1]
     if x < 0 or x >= W:
         raise IndexError(
             f"x coordinate ({x}) if out of bounds for array with width {W}.")
     if y < 0 or y >= H:
         raise IndexError(
             f"y coordinate ({y}) if out of bounds for array with height {H}.")
-        
+
     if img[y, x] == 0:
         raise ValueError(
             f"The pose {pose} corresponds to an occupied cell, cannot raycast from an occupied cell.")
+
+
+def _check_map(img: NDArray) -> NDArray[np.uint8]:
+    if not np.issubdtype(img.dtype, np.integer):
+        raise ValueError(
+            f"Received an array with dtype {img.dtype}, The only supported dtypes are subtypes of np.integral.")
+
+    shape = img.shape
+    if len(shape) < 2 or (len(shape) > 2 and shape[2] > 1):
+        raise ValueError(
+            f"Received an array of shape {shape}. The only supported shapes are (H x W), (H x W x 1)."
+        )
+
+    if len(img.shape) == 3:   # take only the first channel
+        img = img[:, :, 0]
+    if img.dtype != np.uint8:
+        warnings.warn(
+            f"Converting image of dtype {img.dtype} to uint8 for raycasting.")
+        img = img.astype(np.uint8)
+
+    return img
diff --git a/test/test_lidar2D.py b/test/test_lidar2D.py
new file mode 100644
index 0000000..8927aba
--- /dev/null
+++ b/test/test_lidar2D.py
@@ -0,0 +1,194 @@
+import math
+import warnings
+
+import numpy as np
+import pytest
+
+from raycast2D.raycast2D import Lidar2D
+import raycast2D.raycast2D as rc_mod
+
+
+def _scan(
+    img: np.ndarray,
+    pose: tuple[int, int] | tuple[int, int, float],
+    *,
+    num_rays: int = 1000,
+    FOV: int = 360,
+    ray_length: int = 500,
+    only_true_collisions: bool = True,
+):
+    lidar = Lidar2D(num_rays=num_rays, FOV=FOV, ray_length=ray_length)
+    return lidar.scan(pose, image=img, only_true_collisions=only_true_collisions)
+
+
+def test_empty_map_single_ray_hits_border_right() -> None:
+    img = np.full((10, 10), 255, dtype=np.uint8)
+    pose = (5, 5, math.pi / 2)  # yaw=pi/2 -> right (yaw=0 points up)
+
+    rays = _scan(img, pose, num_rays=1, FOV=1, ray_length=100)
+    assert rays.shape == (1, 2)
+
+    x, y = map(int, rays[0])
+    assert (x, y) == (img.shape[1] - 1, pose[1])
+
+
+def test_obstacle_above_origin_collision_is_obstacle_cell() -> None:
+    img = np.full((12, 12), 255, dtype=np.uint8)
+    pose = (6, 9, 0.0)  # up (yaw=0 points up)
+    img[4, 6] = 0  # obstacle directly above along same column
+
+    rays = _scan(img, pose, num_rays=1, FOV=1, ray_length=200)
+    x, y = map(int, rays[0])
+    assert (x, y) == (6, 4)
+
+
+def test_yaw_controls_direction_right_vs_left() -> None:
+    img = np.full((11, 11), 255, dtype=np.uint8)
+    pose = (5, 5)
+
+    rays_right = _scan(img, (pose[0], pose[1], math.pi / 2),
+                       num_rays=1, FOV=1, ray_length=200)
+    xr, yr = rays_right[0]
+    assert (xr, yr) == (img.shape[1] - 1, pose[1])
+
+    rays_left = _scan(img, (pose[0], pose[1], -math.pi / 2),
+                      num_rays=1, FOV=1, ray_length=200)
+    xl, yl = rays_left[0]
+    assert (xl, yl) == (0, pose[1])
+
+
+def test_fov_180_default_yaw_first_left_last_right() -> None:
+    img = np.full((21, 21), 255, dtype=np.uint8)
+    pose = (10, 10)  # default yaw=0 (points up)
+
+    rays = _scan(img, pose, num_rays=2, FOV=180, ray_length=500)
+    (x0, y0) = map(int, rays[0])
+    (x1, y1) = map(int, rays[-1])
+
+    assert (x0, y0) == (0, pose[1])
+    assert (x1, y1) == (img.shape[1] - 1, pose[1])
+
+
+def test_fov_180_yaw_pi_over_2_first_top_last_bottom() -> None:
+    img = np.full((21, 21), 255, dtype=np.uint8)
+    pose = (10, 10, math.pi / 2)
+
+    rays = _scan(img, pose, num_rays=2, FOV=180, ray_length=500)
+    (x0, y0) = map(int, rays[0])
+    (x1, y1) = map(int, rays[-1])
+
+    assert (x0, y0) == (pose[0], 0)
+    assert (x1, y1) == (pose[0], img.shape[0] - 1)
+
+
+def test_pose_on_occupied_cell_raises() -> None:
+    img = np.full((10, 10), 255, dtype=np.uint8)
+    img[3, 4] = 0
+    with pytest.raises(ValueError):
+        _scan(img, (4, 3))
+
+
+def test_rejects_non_integer_dtype() -> None:
+    img = np.zeros((10, 10), dtype=np.float32)
+    img[:] = 1.0
+    with pytest.raises(ValueError):
+        _scan(img, (5, 5))
+
+
+def test_raises_when_image_dtype_is_not_np_integral() -> None:
+    img = np.full((10, 10), 1.0, dtype=np.float64)
+    with pytest.raises(ValueError):
+        _scan(img, (5, 5))
+
+
+def test_accepts_single_channel_3d_image() -> None:
+    img = np.full((10, 10, 1), 255, dtype=np.uint8)
+    rays = _scan(img, (5, 5, 0.0), num_rays=1, FOV=1, ray_length=50)
+    assert rays.shape == (1, 2)
+
+
+def test_only_true_collisions_true_returns_empty_when_no_collision() -> None:
+    img = np.full((30, 30), 255, dtype=np.uint8)
+    pose = (15, 15, 0.0)
+
+    # Small ray_length keeps endpoints inside free space (no border hit, no obstacles)
+    rays_true = _scan(
+        img,
+        pose,
+        num_rays=5,
+        FOV=30,
+        ray_length=3,
+        only_true_collisions=True,
+    )
+    assert rays_true.shape == (0, 2)
+
+    rays_false = _scan(
+        img,
+        pose,
+        num_rays=5,
+        FOV=30,
+        ray_length=3,
+        only_true_collisions=False,
+    )
+    assert rays_false.shape == (5, 2)
+
+
+def test_only_true_collisions_filters_to_obstacle_hits_subset() -> None:
+    img = np.full((30, 30), 255, dtype=np.uint8)
+    pose = (15, 15, math.pi / 2)
+
+    # With num_rays=2 and FOV=180 around yaw=pi/2, angles are exactly -pi/2 (up) and +pi/2 (down)
+    # Put an obstacle only on the "up" ray.
+    img[13, 15] = 0  # (x=15, y=13)
+
+    rays_true = _scan(
+        img,
+        pose,
+        num_rays=2,
+        FOV=180,
+        ray_length=3,
+        only_true_collisions=True,
+    )
+    assert rays_true.shape == (1, 2)
+    assert tuple(map(int, rays_true[0])) == (15, 13)
+
+    rays_false = _scan(
+        img,
+        pose,
+        num_rays=2,
+        FOV=180,
+        ray_length=3,
+        only_true_collisions=False,
+    )
+    assert rays_false.shape == (2, 2)
+
+    pts = {tuple(map(int, p)) for p in rays_false}
+    # One ray hits the obstacle; the other ends in free space (no collision)
+    assert (15, 13) in pts
+    assert (15, 18) in pts
+
+
+def test_uint8_conversion_warning_emitted_only_once() -> None:
+    img = np.full((20, 20), 255, dtype=np.int32)
+    pose = (10, 10, 0.0)
+
+    # Make this test deterministic: clear the module warning registry so the first call warns.
+    registry = getattr(rc_mod, "__warningregistry__", None)
+    if isinstance(registry, dict):
+        registry.clear()
+
+    with warnings.catch_warnings(record=True) as recorded:
+        warnings.simplefilter("default")
+
+        _scan(img, pose, num_rays=1, FOV=1,
+              ray_length=3, only_true_collisions=False)
+        _scan(img, pose, num_rays=1, FOV=1,
+              ray_length=3, only_true_collisions=False)
+
+    conversion_warnings = [
+        w
+        for w in recorded
+        if issubclass(w.category, UserWarning)
+        and "Converting image of dtype" in str(w.message)
+    ]
+    assert len(conversion_warnings) == 1
diff --git a/test/test_raycast2D.py b/test/test_raycast2D.py
index f7a89e3..163ba2b 100644
--- a/test/test_raycast2D.py
+++ b/test/test_raycast2D.py
@@ -10,7 +10,7 @@
 
 def test_empty_map_single_ray_hits_border_right() -> None:
     img = np.full((10, 10), 255, dtype=np.uint8)
-    pose = (5, 5, 0.0)  # yaw=0 -> right
+    pose = (5, 5, math.pi / 2)  # yaw=pi/2 -> right (yaw=0 points up)
 
     rays = cast(img, pose, num_rays=1, FOV=1, ray_length=100)
     assert rays.shape == (1, 2)
@@ -21,7 +21,7 @@ def test_empty_map_single_ray_hits_border_right() -> None:
 
 def test_obstacle_above_origin_collision_is_obstacle_cell() -> None:
     img = np.full((12, 12), 255, dtype=np.uint8)
-    pose = (6, 9, -math.pi / 2)  # up
+    pose = (6, 9, 0.0)  # up (yaw=0 points up)
     img[4, 6] = 0  # obstacle directly above along same column
 
     rays = cast(img, pose, num_rays=1, FOV=1, ray_length=200)
@@ -33,12 +33,12 @@ def test_yaw_controls_direction_right_vs_left() -> None:
     img = np.full((11, 11), 255, dtype=np.uint8)
     pose = (5, 5)
 
-    rays_right = cast(img, (pose[0], pose[1], 0.0),
+    rays_right = cast(img, (pose[0], pose[1], math.pi / 2),
                       num_rays=1, FOV=1, ray_length=200)
     xr, yr = rays_right[0]
     assert (xr, yr) == (img.shape[1] - 1, pose[1])
 
-    rays_left = cast(img, (pose[0], pose[1], -math.pi),
+    rays_left = cast(img, (pose[0], pose[1], -math.pi / 2),
                      num_rays=1, FOV=1, ray_length=200)
     xl, yl = rays_left[0]
     assert (xl, yl) == (0, pose[1])
@@ -46,7 +46,7 @@ def test_yaw_controls_direction_right_vs_left() -> None:
 
 def test_fov_180_default_yaw_first_left_last_right() -> None:
     img = np.full((21, 21), 255, dtype=np.uint8)
-    pose = (10, 10)  # default yaw=-pi/2
+    pose = (10, 10)  # default yaw=0 (points up)
 
     rays = cast(img, pose, num_rays=2, FOV=180, ray_length=500)
     (x0, y0) = map(int, rays[0])
@@ -56,9 +56,9 @@ def test_fov_180_default_yaw_first_left_last_right() -> None:
     assert (x1, y1) == (img.shape[1] - 1, pose[1])
 
 
-def test_fov_180_yaw0_first_top_last_bottom() -> None:
+def test_fov_180_yaw_pi_over_2_first_top_last_bottom() -> None:
     img = np.full((21, 21), 255, dtype=np.uint8)
-    pose = (10, 10, 0.0)
+    pose = (10, 10, math.pi / 2)
 
     rays = cast(img, pose, num_rays=2, FOV=180, ray_length=500)
     (x0, y0) = map(int, rays[0])
@@ -122,9 +122,9 @@ def test_only_true_collisions_true_returns_empty_when_no_collision() -> None:
 
 def test_only_true_collisions_filters_to_obstacle_hits_subset() -> None:
     img = np.full((30, 30), 255, dtype=np.uint8)
-    pose = (15, 15, 0.0)
+    pose = (15, 15, math.pi / 2)
 
-    # With num_rays=2 and FOV=180 around yaw=0, angles are exactly -pi/2 (up) and +pi/2 (down)
+    # With num_rays=2 and FOV=180 around yaw=pi/2, angles are exactly -pi/2 (up) and +pi/2 (down)
     # Put an obstacle only on the "up" ray.
     img[13, 15] = 0  # (x=15, y=13)