🧭 Concept LiDAR & point clouds Files & formats LiDAR & Point Clouds · step 2 of 14 Updated

LAS, LAZ and COPC Explained: How Point Clouds Are Stored

Problem statement

A point cloud file is a header plus a long run of fixed-width records. That layout makes it fast to read from start to finish and slow to answer the question everyone actually asks: give me the points in this box.

The formats exist on a ladder that trades compression and indexing against simplicity:

  • LAS — the uncompressed base format. Simple, large, no spatial index.
  • LAZ — LAS compressed in chunks. Typically 5–10× smaller, still no spatial index.
  • COPC — LAZ with an octree index inside it, so a bounding-box query at a chosen resolution is a range request.
  • EPT — an octree stored as a directory of many small LAZ files plus a JSON hierarchy.

Measured on a real 3DEP survey: 12,968,770 points over 877 × 876 m arrived as 287 LAZ files totalling 82.8 MB — 6.4 bytes per point on disk, against 34 bytes per point for the same records uncompressed.

Quick answer 

import laspy

las = laspy.read("tile.laz")
h = las.header
print(f"LAS {h.version}, point format {h.point_format.id}")
print(f"{h.point_count:,} points")
print(f"scales  {h.scales}")
print(f"offsets {h.offsets}")
print(f"raw X range {las.X.min()} .. {las.X.max()}   (int32)")
print(f"scaled x    {las.x.min():.2f} .. {las.x.max():.2f}")
LAS 1.2, point format 1
180,756 points
scales  [0.01 0.01 0.01]
offsets [-1.73345e+07  8.07826e+06  4.86000e+02]
raw X range -802212 .. -640954   (int32)
scaled x    -17342522.12 .. -17340909.54

las.X and las.x are different attributes: the raw 32-bit integer and the scaled coordinate. Using the wrong one is off by a factor of 100 here, and by 1,000 in files with millimetre scaling.

A LAS file as a header carrying scales and offsets, followed by fixed-width point records storing coordinates as scaled integers.
Coordinates are integers. The header holds the scale and offset that turn them into positions.

Step-by-step solution

1. Understand why coordinates are integers

LAS stores X, Y and Z as signed 32-bit integers, with a scale and offset per axis in the header:

position = raw_integer * scale + offset

That is a deliberate design. Integers of a fixed scale give uniform precision across the file — every coordinate is exact to 1 cm here — whereas floating point gives relative precision, so a coordinate of 17 million metres would lose resolution exactly where a projected CRS puts you.

The cost is a fixed range. With a 0.01 scale, a signed 32-bit integer spans about ±21 million, which is why the offset exists: it moves the origin near the data.

2. Know what the point format determines

The point format ID fixes which attributes exist. Format 1 adds GPS time to the basic format 0; formats 2 and 3 add RGB; formats 6 and above (LAS 1.4) extend the return numbering and the classification range.

print([d.name for d in las.point_format.dimensions])
['X', 'Y', 'Z', 'intensity', 'return_number', 'number_of_returns',
 'scan_direction_flag', 'edge_of_flight_line', 'classification',
 'synthetic', 'key_point', 'withheld', 'scan_angle_rank', 'user_data',
 'point_source_id', 'gps_time', 'OriginId']

An older format cannot hold more than 5 returns per pulse or classification codes above 31. If you need those, the file must be LAS 1.4 with point format 6 or higher.

3. Know what LAZ does and does not add

LAZ compresses the point records in independent chunks, typically 50,000 points each. Compression is lossless and typically gets a factor of 5–10 — here 82.8 MB for records that would occupy about 440 MB uncompressed.

What LAZ does not add is a spatial index. Finding the points in a box still means decompressing every chunk and testing. The chunking helps only in that a reader can skip whole chunks whose declared bounds miss the query, if the writer recorded them.

4. Use COPC or EPT when you need spatial queries

COPC stores the points in octree order inside a single LAZ file, with a hierarchy that maps each octree node to a byte range. A viewer can then fetch exactly the nodes overlapping the view, at exactly the depth the zoom level needs, using HTTP range requests.

EPT stores the same idea as a directory: ept.json, a hierarchy of JSON files, and one LAZ per node. The survey measured here is EPT, which is why it arrived as 287 separate files.

depth 0:      1 node,      19,873 points
depth 3:     31 nodes,   1,857,900 points
depth 5:    317 nodes,  44,739,745 points
depth 7:  4,720 nodes, 281,081,360 points

Each depth is a complete, coarser sample of the whole cloud — so a low-resolution overview costs one small node rather than the whole file.

5. Check the CRS, and do not assume it is metric

The survey measured here is stored in EPSG:3857, Web Mercator. Its extent reads as 1,613 × 1,613 m and the ground truth is 877 × 876 m — a factor of 1.92 at 58.6° north.

Every metric quantity computed in those units is wrong. See laspy returns the wrong coordinates.

LAS, LAZ, COPC and EPT compared on compression, spatial indexing and whether a bounding-box query is cheap.
LAZ compresses; COPC and EPT index. Compression without an index still means reading everything.

Code examples

Example 1 — a header report before you read any points 

import laspy
import numpy as np


def header_report(path):
    """Everything the header knows, without decompressing the points."""
    with laspy.open(path) as reader:
        h = reader.header
        print(f"  LAS {h.version}, point format {h.point_format.id}")
        print(f"  {h.point_count:,} points")
        print(f"  scales  {h.scales}")
        print(f"  offsets {h.offsets}")
        print(f"  bounds  x {h.mins[0]:.2f}..{h.maxs[0]:.2f}")
        print(f"          y {h.mins[1]:.2f}..{h.maxs[1]:.2f}")
        print(f"          z {h.mins[2]:.2f}..{h.maxs[2]:.2f}")

        crs = h.parse_crs()
        print(f"  CRS     {crs.to_string() if crs else 'NOT SET'}")
        if crs and not crs.is_projected:
            print("  ! geographic CRS — reproject before any metric computation")

        extent = (h.maxs[0] - h.mins[0]) * (h.maxs[1] - h.mins[1])
        print(f"  nominal density {h.point_count / extent:.1f} per square unit"
              f" (units of the file CRS)")
        return h

Opening with laspy.open reads the header only. On an 80 MB file that is instant, against seconds to decompress every point — worth it when you are surveying a directory of tiles.

Example 2 — reading in chunks, so a large file fits 

import laspy
import numpy as np


def chunked_stats(path, chunk_size=1_000_000, ground_class=2):
    """Stream a large LAZ without materialising it."""
    total = 0
    z_min, z_max = np.inf, -np.inf
    class_counts = {}

    with laspy.open(path) as reader:
        for points in reader.chunk_iterator(chunk_size):
            total += len(points)
            z = np.asarray(points.z)
            z_min = min(z_min, z.min())
            z_max = max(z_max, z.max())
            codes, counts = np.unique(points.classification, return_counts=True)
            for code, count in zip(codes.tolist(), counts.tolist()):
                class_counts[code] = class_counts.get(code, 0) + count

    print(f"  {total:,} points, z {z_min:.2f}..{z_max:.2f}")
    for code, count in sorted(class_counts.items(), key=lambda kv: -kv[1]):
        print(f"    class {code:2d}: {count:12,}  {count / total:6.2%}")
    return class_counts

A 13-million-point cloud is 311 MB as float64 coordinates alone. Streaming keeps memory flat, and for anything that reduces — histograms, grids, counts — costs nothing.

Example 3 — merging a directory of tiles safely 

import glob
import laspy
import numpy as np


def merge_tiles(pattern, max_points=None):
    """Concatenate tiles, checking that they agree on scale and CRS."""
    files = sorted(glob.glob(pattern))
    if not files:
        raise FileNotFoundError(pattern)

    with laspy.open(files[0]) as reader:
        reference = reader.header
        scales, crs = reference.scales, reference.parse_crs()

    xs, ys, zs, cls, ret, nret = [], [], [], [], [], []
    total = 0
    for path in files:
        las = laspy.read(path)
        if not np.allclose(las.header.scales, scales):
            raise ValueError(f"{path}: scales {las.header.scales} "
                             f"differ from {scales}")
        if las.header.parse_crs() != crs:
            raise ValueError(f"{path}: CRS differs from the first tile")
        xs.append(las.x); ys.append(las.y); zs.append(las.z)
        cls.append(np.asarray(las.classification))
        ret.append(np.asarray(las.return_number))
        nret.append(np.asarray(las.number_of_returns))
        total += len(las.points)
        if max_points and total >= max_points:
            break

    print(f"  merged {len(xs)} tiles, {total:,} points")
    return (np.concatenate(xs), np.concatenate(ys), np.concatenate(zs),
            np.concatenate(cls), np.concatenate(ret), np.concatenate(nret))
  merged 287 tiles, 12,968,770 points

Checking the scales matters more than it looks. Tiles written with different scales still merge without error in x (scaled) space, and produce a cloud that cannot be written back to a single LAS file without re-scaling.

Explanation

Why fixed-point integers instead of floats

Floating point has relative precision: a float32 has about 7 significant digits, so at a UTM easting of 500,000 m it resolves about 3 cm, and at a Web Mercator easting of 17,000,000 m about 1 m.

Scaled integers have absolute precision. Every coordinate in the file measured here is exact to 1 cm regardless of position, because the offset moves the origin next to the data and the integer counts centimetres from there.

The trade is a bounded range and a mandatory header lookup. A reader that ignores the scale and offset gets numbers that are wrong by a factor of 100 and look plausible.

Why LAZ does not make queries fast

Compression reduces bytes; it does not reorganise. In a LAZ file the points are still in acquisition order, so the points inside any bounding box are scattered through every chunk.

The consequence is that "read a small area from a big file" costs almost as much as reading the whole file. That is exactly the operation a web viewer, a tiling job or an interactive tool needs, which is why the indexed formats exist.

Why an octree is the right index

Point clouds need two things from an index: spatial selection and level of detail. An octree gives both from one structure.

Each node holds a sample of the points in its cube, and its children hold finer samples of the eight sub-cubes. Reading to depth 3 gives a coarse view of everything; descending only where the viewport is gives detail where it is needed.

That is why the depth counts above rise to a peak and then fall: the deepest level holds only the leftover points in the densest regions, not a uniform layer.

Why merging tiles is not always safe

Adjacent tiles from one survey normally share a CRS, scale and offset. Tiles from different deliveries, or reprocessed at different times, may not.

The failure is silent. Scaled coordinates merge fine, and only the write-back fails — or worse, succeeds with a re-scaled precision that quietly loses a decimal place. Assert on scales and parse_crs() at merge time.

EPT octree depths holding 19,873 points at depth 0 rising to 281 million at depth 7 and falling at depth 8.
The count falls at the deepest level because it holds only the leftover points in dense regions.

Edge cases or notes

  • las.X is the raw integer; las.x is scaled. Off by 100 or 1,000.
  • The offset is not the minimum. It is an arbitrary origin near the data.
  • Point format limits the attributes. Formats below 6 cap returns at 5 and classes at 31.
  • LAZ has no spatial index. Use COPC or EPT for bounding-box queries.
  • CRS may be missing entirely. parse_crs() returning None is common in older files.
  • Check the CRS is projected before computing anything metric.
  • Stream large files with chunk_iterator; 13 M points is 311 MB of coordinates alone.
  • Assert scale and CRS agreement when merging tiles.

FAQ

What is the difference between LAS and LAZ?

LAZ is losslessly compressed LAS, typically 5–10× smaller. It adds no spatial index, so bounding-box queries still cost a full read.

What is COPC?

Cloud-Optimised Point Cloud: a LAZ file whose points are stored in octree order with a byte-range index, so a viewer can fetch just the region and resolution it needs over HTTP.

Why are point cloud coordinates stored as integers?

For uniform absolute precision. raw * scale + offset gives every point the same resolution regardless of position, which floating point cannot at large coordinate values.

Why is las.X different from las.x?

X is the raw 32-bit integer as stored; x applies the header's scale and offset. Confusing them is off by a factor of 100 or 1,000.

How do I read a huge LAZ file without running out of memory?

Use laspy.open(...).chunk_iterator(n) and reduce as you go. A 13-million-point cloud is 311 MB of coordinates alone as float64.

Can I merge LAS tiles by concatenating them?

Only after checking that they share a CRS, scale and offset. Mismatched scales merge silently in scaled space and fail or lose precision on write-back.

What is EPT?

Entwine Point Tile: an octree stored as a directory of small LAZ files plus a JSON hierarchy. The 3DEP survey used here arrived as 287 such files.