Introduction
This article uses Depth Anything 3 (DA3) to perform monocular depth estimation on a single RGB GeoTIFF. No additional training is performed; we use a publicly available pretrained model.
Depth estimated from a single satellite image can be used to generate building models with geometry comparable to LOD2. For an overview of monocular depth estimation and our work on generating 3D building data from satellite imagery, see this explanatory page.
As the first step in that workflow, this article covers environment setup on Ubuntu with an NVIDIA GPU, model downloading, depth inference, GeoTIFF export, and inspection in QGIS. Building polygons, a DEM, and camera exterior orientation data are not used as inputs.
For the satellite imagery example, we used an archived WorldView-2 image with a resolution of 50 cm. Because lossy compression noise made building boundaries unclear in the original image, we first applied 2× super-resolution and image enhancement.
The resulting GeoTIFF is 2048×1955 pixels. When we performed monocular depth estimation on the entire image, the predicted boundaries appeared blurred, possibly because the whole image was downscaled internally. We therefore tried splitting the image into 300-pixel tiles, running inference on each tile, and blending the overlapping predictions.
This approach was inspired by our experience with Depth Anything 2 a little over a year earlier, when depth estimation for the same area produced a map that appeared tilted across the entire scene.
For our image, tiling produced visually sharper building boundaries, as it had with Depth Anything 2. However, this is a visual comparison, not a validation of building heights or roof slopes against measurements.
The output is relative depth. It does not directly provide building heights in meters, a DSM, or a DEM. Generating 3D building models and placing them on a DEM are separate postprocessing steps.
1. Model and Processing Pipeline
We use depth-anything/DA3MONO-LARGE. Its official model card describes it as a model for monocular relative depth estimation. It is a different model from DA3-LARGE and DA3METRIC-LARGE.
The workflow is: read the RGB GeoTIFF → run DA3 inference → resize the depth prediction to the original pixel grid → save as GeoTIFF. The tiled version splits the image before inference and blends the depth tiles before saving.
We use the official API for DA3 itself. GeoTIFF input/output and tile blending are implemented in the additional scripts provided in this article. This is not an official DA3 geospatial tiled inference feature.
2. Environment
The environment used in this article is listed below. All commands are intended for Bash on Ubuntu.
| Item | Configuration |
|---|---|
| OS | Ubuntu 24.04 |
| GPU | NVIDIA RTX PRO 6000 Blackwell Max-Q |
| Python | 3.10 / dedicated conda environment |
| PyTorch | 2.9.1, CUDA 12.8 wheel |
| torchvision | 0.24.1 |
| xFormers | 0.0.33.post2 |
| NumPy | 1.26.4 |
| OpenCV | opencv-python 4.10.0.84 |
| MoviePy | 1.0.3 |
| Model | DA3MONO-LARGE |
| GeoTIFF input/output | Rasterio |
The memory capacity of this GPU is not a requirement for using the workflow. GPU memory requirements depend on input resolution and the environment. This article does not evaluate minimum memory requirements or processing speed on other GPUs.
For the PyTorch version combination, see the official previous versions page. For xFormers compatibility, see the 0.0.33.post2 release notes. This procedure pins a combination of major packages rather than upgrading everything to the latest version.
The “CUDA Version” shown by nvidia-smi, the Toolkit version shown by nvcc --version, and the CUDA version used by PyTorch, shown by torch.version.cuda, have different meanings. This procedure uses the CUDA 12.8 PyTorch wheel and does not replace the system CUDA Toolkit. Install a suitable NVIDIA GPU driver beforehand.
3. Create the Working Directory and Conda Environment
We assume that Miniconda or Anaconda is already available.
conda create -n da3_bw python=3.10 pip -y
conda activate da3_bw
mkdir -p ~/da3_geotiff/input ~/da3_geotiff/output ~/da3_geotiff/logs
cd ~/da3_geotiff
nvidia-smi
git --version
tmux -V
If git or tmux is unavailable, install it through Ubuntu's package manager.
sudo apt-get update
sudo apt-get install -y git tmux
Save all files listed in “Appendix: Scripts” at the end of this article directly under ~/da3_geotiff/. Copy the code from each expandable section and save it using the filename in its heading. No separate custom ZIP download or private repository is required.
| File | Purpose |
|---|---|
constraints.txt |
Pin major package versions |
setup.sh |
Install official code and dependencies |
check_env.py |
Check CUDA and imports |
download_model.py |
Download official model weights |
infer_geotiff.py |
Run inference on the whole image |
da3_tile_io.py |
Handle GeoTIFF input/output for tiled inference |
infer_geotiff_tiled.py |
Run tiled inference and blend predictions |
4. Install DA3 and Its Dependencies
conda activate da3_bw
cd ~/da3_geotiff
bash setup.sh
setup.sh performs the following steps:
- Install the CUDA 12.8 versions of PyTorch, torchvision, and xFormers.
- Clone the official repository into
vendor/Depth-Anything-3. - Install the basic DA3 dependencies, Rasterio, and other required packages.
- Run
pip check, check imports, and test matrix multiplication on CUDA. - Save the Git commit ID and
pip freezeoutput inlogs/.
The official code is obtained with the following command:
# This is executed by setup.sh; you do not need to run it again manually.
git clone https://github.com/ByteDance-Seed/Depth-Anything-3.git vendor/Depth-Anything-3
Only the basic dependencies are installed; optional components for Gradio and 3D Gaussian output are not installed. Depth inference uses infer_gs=False.
To check the environment again, run:
python -m pip check
python check_env.py
Check that the output includes CUDA runtime: 12.8, the GPU name, and Imports and CUDA matrix multiplication: OK. These checks alone do not establish that inference with the model weights succeeds.
Although major packages are pinned, this is not a fully locked environment with every dependency and repository commit fixed in advance. If future changes cause dependency conflicts, inspect the error log and downloaded commit before removing the constraints and upgrading packages.
5. Download the DA3MONO-LARGE Weights
python download_model.py
The model is downloaded from the official Hugging Face repository, depth-anything/DA3MONO-LARGE, into models/DA3MONO-LARGE/. The script resolves a revision, downloads that specific revision, and records its ID in download_revision.txt.
cat models/DA3MONO-LARGE/download_revision.txt
cat logs/da3_git_commit.txt
These files record the versions retrieved at download time. They do not guarantee that a new download on another day will resolve to the same revision. Retain the model directory and environment records for reproducibility.
Inference loads the locally saved weights. The scripts do not upload input imagery to Hugging Face.
6. Prepare the Input GeoTIFF
The examples use input/area_rgb.tif. Place your own image at that path, or change the input path in the commands below.
The input requirements are:
- An RGB source image with an 8-bit
uint8pixel type - A GeoTIFF with a CRS and affine transform
- Bands 1, 2, and 3 used as R, G, and B by default
- No previously estimated depth or height image as input
For 16-bit satellite imagery, first convert it to 8-bit RGB using a consistent method across the entire image. Avoid a simple astype(uint8) conversion or independent contrast adjustments for each tile. The script stops on 16-bit imagery rather than converting it automatically.
Example input inspection:
python - <<'PY'
import rasterio
with rasterio.open('input/area_rgb.tif') as src:
print('size:', src.width, src.height)
print('bands:', src.count)
print('dtype:', src.dtypes)
print('CRS:', src.crs)
print('transform:', src.transform)
print('color interpretation:', src.colorinterp)
PY
If the RGB band order differs, add --bands R_number G_number B_number when running inference. For example, use --bands 3 2 1 if R, G, and B correspond to bands 3, 2, and 1.
An additional band in a four-band image is not passed to the model. However, alpha and NoData masks recognized by Rasterio are used to identify valid pixels. Missing pixels are filled with the RGB median during inference and restored to NoData in the output. This does not eliminate their influence on predictions near missing regions.
7. Start with Whole-Image Inference
Use tmux for long-running processing.
tmux new -s da3
Run the following commands inside the session:
conda activate da3_bw
cd ~/da3_geotiff
mkdir -p output logs
set -o pipefail
python -u infer_geotiff.py \
--input input/area_rgb.tif \
--output-prefix output/area_da3_full_r1008 \
--process-res 1008 \
2>&1 | tee logs/area_da3_full_r1008.log
--process-res 1008 specifies the approximate longest side of the image passed to the model. It does not split the original image into 1008-pixel tiles. Preprocessing resizes the image without cropping and adjusts the resolution to the model's patch size.
The output depth is resized to the original image dimensions using bilinear interpolation before saving. Therefore, matching the input and output pixel dimensions does not mean that details lost during downscaled inference have been recovered.
Press Ctrl+B, then D, to detach from tmux. Reattach with:
tmux attach -t da3
8. Use 300-Pixel Tiles to Inspect Building Details
Instead of downscaling the whole image, run inference on small regions of the original image. Here we use 300-pixel crops, 100-pixel overlaps, and a model input size of 518 pixels.
python -u infer_geotiff_tiled.py \
--input input/area_rgb.tif \
--output-prefix output/area_da3_t300_o100_r518_shift \
--tile-size 300 \
--overlap 100 \
--process-res 518 \
--alignment shift \
2>&1 | tee logs/area_da3_t300_o100_r518_shift.log
| Argument | Meaning |
|---|---|
--tile-size 300 |
Crop size in the original image, normally 300×300 pixels |
--overlap 100 |
Overlap between adjacent tiles; the normal stride is 200 pixels |
--process-res 518 |
Input size used when passing a 300-pixel RGB tile to the model |
--alignment shift |
Correct only the offset using depth differences in overlapping regions |
--output-prefix |
Common filename prefix for the output files |
518 is a multiple of 14. This wrapper requires process-res to be a multiple of 14. Upscaling a 300-pixel tile to 518 pixels is a model input preparation step; it does not improve the actual ground resolution of the imagery.
For example, a 300-pixel tile covers 75×75 m in an image with a resolution of 0.25 m/pixel. At 0.5 m/pixel, it covers 150×150 m. The proportion occupied by buildings and the surrounding context therefore differ even with the same 300-pixel tile size. Equivalent accuracy is not guaranteed for other regions or image resolutions.
At image boundaries, crop positions are shifted back to keep the tiles inside the image, so the final overlap may exceed 100 pixels. If the image width or height is less than 300 pixels, the image is padded using edge pixel values. A 2048×1955 image produces 10×10 tiles, or 100 tiles in total.
The model is loaded once at the start, and tiles are processed sequentially. Each prediction is resized back to the original 300-pixel tile grid and blended into the mosaic.
Blending Modes: shift, affine, and none
| Mode | Correction | Characteristics |
|---|---|---|
shift |
Offset only | Does not rescale depth differences within each tile |
affine |
Scale and offset | Aligns scales in overlapping regions, but also changes the magnitude of local relief |
none |
No value correction | Takes a weighted average of the original predictions; differences between tiles may remain |
In this workflow, we explicitly use shift to avoid the effects of scale correction on local height differences during comparison. The supplied code defaults to affine, so do not omit this argument.
For shift alignment, let be the current tile's depth and the depth in the overlapping region of the existing mosaic. For valid overlapping pixels, apply:
For the first tile, or a tile without valid overlap, the offset is set to zero. The corrected tiles are then blended with a weighted average that assigns greater weights toward each tile's center.
“Alignment” here does not geometrically transform the image. It aligns the reference levels of depth values predicted at the same pixel locations.
The tiles are not independently normalized to 0–1 using their own minimum and maximum values before blending. Nevertheless, shift alignment cannot resolve true scale differences. Scan-order-dependent error accumulation and blurring caused by averaging overlapping regions may remain.
9. Run Inference on Another Area from Scratch
Change the input path and output name, then run the same command. Depth results or tiles from the previous area are not reused.
python -u infer_geotiff_tiled.py \
--input input/another_area_rgb.tif \
--output-prefix output/another_area_da3_t300_o100_r518_shift \
--tile-size 300 --overlap 100 --process-res 518 \
--alignment shift \
2>&1 | tee logs/another_area_da3_t300_o100_r518_shift.log
The scripts prevent overwriting existing outputs. When rerunning inference on the same image, specify an unused --output-prefix. The DA3 scripts in this article do not support resuming an interrupted run.
10. Understanding the Output Files
| Filename suffix | Contents | Main use |
|---|---|---|
_depth.tif |
Float32 relative depth; larger values indicate greater distance | Numerical analysis and downstream processing |
_near01.tif |
Display GeoTIFF normalized to 0–1, with nearer regions set toward 1 | Visualization in QGIS |
_near.png |
8-bit PNG using the same display representation | Image comparisons and article illustrations |
_metadata.json |
Inference settings, GPU, image dimensions, tile corrections, and other information | Checking execution conditions |
The GeoTIFF outputs retain the input CRS, affine transform, width, and height. The output NoData value is -9999. Invalid pixel locations are preserved, but arbitrary TIFF tags and the original NoData value are not all copied unchanged.
For tiled inference, _depth.tif contains relative depth after blending corrections, not the uncorrected predictions for individual tiles. Shift or affine corrections can produce negative values. The script does not save unblended tile predictions as separate files.
The near01 Display Is Not Building Height
Let be the 1st and 99th percentiles of the valid depth values across the blended mosaic. The display image is computed as:
Nearer regions appear white and farther regions black. This is contrast adjustment and a linear reversal of the depth direction, not inverse depth or a conversion to height in meters.
Because the display range is determined independently for each image, the same brightness in near01 outputs from different regions does not represent the same depth or height. Use _depth.tif for numerical processing.
11. Inspect the Results in QGIS
- Load the original RGB GeoTIFF and
_near01.tifinto QGIS. - Display
_near01.tifusing “Singleband gray.” - Set the minimum to 0, maximum to 1, and color ramp from black to white.
- At the same location and scale, inspect building boundaries, roof ridges, roof surfaces, and their relationship to roads.
When displaying _depth.tif, stretch it according to its actual value range. A fixed 0–1 range may cause many pixels to appear solid white or black. A white image alone does not establish that inference failed.
Example checks for coordinates and output values:
python - <<'PY'
import numpy as np
import rasterio
source='input/area_rgb.tif'
result='output/area_da3_t300_o100_r518_shift_depth.tif'
with rasterio.open(source) as src, rasterio.open(result) as dst:
assert src.crs == dst.crs
assert src.transform == dst.transform
assert (src.width, src.height) == (dst.width, dst.height)
source_valid=np.all(src.read_masks([1,2,3]) > 0, axis=0)
depth=dst.read(1, masked=True)
assert np.array_equal(source_valid, ~np.ma.getmaskarray(depth))
values=depth.compressed()
assert values.size and np.isfinite(values).all()
print('CRS, pixel grid, and valid pixel locations: OK')
print('NoData:', dst.nodata)
print('min / median / max:', values.min(), np.median(values), values.max())
PY
If you change the RGB bands, update [1,2,3] in the checking code accordingly. These checks validate data positioning and format; they do not evaluate depth accuracy.
Example Inference Results
The image below shows the input and the monocular depth estimation result from Depth Anything 3. The input is an archived WorldView-2 image processed with 2× super-resolution.
Input image and Depth Anything 3 monocular depth estimation result:
Left: input image after 2× super-resolution. Right: estimated relative depth displayed with nearer regions brighter.
Differences between buildings and the surrounding ground, as well as relief on some roof surfaces, appear as variations in brightness. However, broad brightness variations are also visible across the image. These may include terrain effects and errors from estimation or tile blending; this image alone cannot distinguish their causes.
Roof brightness should therefore not be interpreted directly as building height. Relative differences from the surrounding ground and agreement with measured values must be checked separately.
12. Common Problems
| Symptom | Checks and actions |
|---|---|
CUDA GPU required / CUDA is not detected |
Check nvidia-smi, python check_env.py, and the active conda environment |
| CUDA out of memory | Lower process-res for whole-image inference or use tiled inference. For tiled inference, compare 518 with another multiple of 14, such as 392. |
xFormers version mismatch |
Check the PyTorch/xFormers combination. Do not upgrade only one of them in the environment. |
Model directory missing |
Complete python download_model.py. If the model is stored elsewhere, specify --model. |
Input must be RGB uint8 |
Check that the input is the source RGB image and that band selection and data type are correct |
| Output appears entirely white or black | Check the display range. Start by displaying near01 over 0–1 and inspecting depth minimum and maximum values. |
| Depth discontinuities at tile boundaries | Compare overlap widths and blending modes. Distinguish actual height differences from blending errors. |
| Buildings appear too high or too low overall | Check whether relative depth is being interpreted as height in meters |
Output exists |
Run from scratch with an unused output prefix |
This implementation tiles model inputs, but it loads the original RGB image and mosaic arrays into CPU memory. It is not a streaming implementation that processes arbitrarily large GeoTIFFs with constant memory usage. Crop an AOI when working with large-area imagery.
13. Separate the Evaluation Questions
Visually sharp building boundaries, correct local depth geometry, and building heights that agree with measurements are separate evaluation questions.
For comparisons, use the same source image, crop extent, and display scale when examining whole-image versus tiled inference and different blending modes. Comparisons with DA2 or other models also require attention to internal input resolution and the representation of output values.
GeoTIFF export in this workflow preserves image georeferencing. It does not calibrate depth to ground-relative height or elevation. To proceed to a DSM or 3D building models, depth-to-height relationships, ground treatment, and validation against measured values must be addressed separately.
References
- Depth Anything 3 official repository
- DA3MONO-LARGE model card
- PyTorch: Previous Versions
- xFormers 0.0.33.post2: wheels for PyTorch 2.9.1
Appendix: Scripts
The following scripts are the GeoTIFF wrappers used in this article. Save each code block directly in the working directory using the indicated filename. The official model implementation and weights are not included here; obtain them from the official sources using the steps above.
constraints.txt
Full contents of constraints.txt
torch==2.9.1
torchvision==0.24.1
xformers==0.0.33.post2
numpy==1.26.4
opencv-python==4.10.0.84
moviepy==1.0.3
setup.sh
Full contents of setup.sh
#!/usr/bin/env bash
set -euo pipefail
cd "$(dirname "$0")"
if [[ "${CONDA_DEFAULT_ENV:-}" != "da3_bw" ]]; then
echo 'Activate the dedicated environment first: conda activate da3_bw' >&2
exit 1
fi
mkdir -p vendor models output logs
python -m pip install --upgrade pip
python -m pip install -c constraints.txt torch torchvision xformers \
--index-url https://download.pytorch.org/whl/cu128
if [[ ! -d vendor/Depth-Anything-3 ]]; then
git clone https://github.com/ByteDance-Seed/Depth-Anything-3.git vendor/Depth-Anything-3
fi
python -m pip install -c constraints.txt -e ./vendor/Depth-Anything-3
python -m pip install -c constraints.txt rasterio addict
python -m pip check
git -C vendor/Depth-Anything-3 rev-parse HEAD > logs/da3_git_commit.txt
python -m pip freeze > logs/environment_freeze.txt
python check_env.py
check_env.py
Full contents of check_env.py
import torch
import torchvision
import xformers
import rasterio
from depth_anything_3.api import DepthAnything3
print('torch:', torch.__version__)
print('torchvision:', torchvision.__version__)
print('xformers:', xformers.__version__)
print('CUDA runtime:', torch.version.cuda)
print('rasterio:', rasterio.__version__)
assert torch.cuda.is_available(), 'CUDA GPU is unavailable'
print('GPU:', torch.cuda.get_device_name(0))
print('Capability:', torch.cuda.get_device_capability(0))
x = torch.randn(512, 512, device='cuda')
y = x @ x
torch.cuda.synchronize()
assert torch.isfinite(y).all()
print('Imports and CUDA matrix multiplication: OK (model inference not yet tested)')
download_model.py
Full contents of download_model.py
from pathlib import Path
from huggingface_hub import HfApi, snapshot_download
root = Path(__file__).resolve().parent
repo = 'depth-anything/DA3MONO-LARGE'
info = HfApi().model_info(repo)
path = snapshot_download(
repo_id=repo, revision=info.sha,
local_dir=str(root / 'models' / 'DA3MONO-LARGE'),
allow_patterns=['*.json', '*.safetensors', '*.md', 'LICENSE*'],
)
(root / 'models' / 'DA3MONO-LARGE' / 'download_revision.txt').write_text(
f'{repo}\n{info.sha}\n', encoding='utf-8')
print('Model downloaded:', path)
infer_geotiff.py
Full contents of infer_geotiff.py
"""Single-image DA3MONO inference. Relative depth, never metric height.
Whole-image inference avoids independent tile scales. The API resizes, without
cropping; the depth is bilinearly restored to the source pixel grid. Original
CRS/affine and the intersection of RGB valid masks are retained. Unknown extra
bands are ignored; rasterio-recognized alpha masks are honored.
"""
from pathlib import Path
import argparse
import json
import time
import numpy as np
import rasterio
from PIL import Image
NODATA = -9999.0
def read_rgb(path, bands):
with rasterio.open(path) as src:
if min(bands) < 1 or max(bands) > src.count:
raise ValueError(f'Invalid RGB bands {bands}; source has {src.count} bands')
rgb = src.read(bands)
if rgb.dtype != np.uint8:
raise ValueError('Input must be RGB uint8. Convert 16-bit/float data explicitly first.')
if src.crs is None:
raise ValueError('Input has no CRS; supply a georeferenced RGB GeoTIFF')
valid = np.all(src.read_masks(bands) > 0, axis=0)
if not valid.any():
raise ValueError('Input has no valid RGB pixels')
profile = dict(driver='GTiff', width=src.width, height=src.height,
count=1, dtype='float32', crs=src.crs, transform=src.transform,
nodata=NODATA, compress='deflate', predictor=3, tiled=True,
blockxsize=256, blockysize=256, BIGTIFF='IF_SAFER')
rgb = np.moveaxis(rgb, 0, -1).copy()
# Fill masked pixels only for inference; their outputs remain NoData.
if not valid.all():
rgb[~valid] = np.median(rgb[valid], axis=0).astype(np.uint8)
print('Warning: masked areas can affect nearby predictions.', flush=True)
return rgb, valid, profile
def restore_depth(depth, height, width):
depth = np.asarray(depth, dtype=np.float32)
if depth.ndim != 2 or not np.isfinite(depth).all():
raise ValueError('Model depth must be a finite 2-D array')
if depth.shape != (height, width):
depth = np.asarray(Image.fromarray(depth).resize(
(width, height), resample=Image.Resampling.BILINEAR), dtype=np.float32)
return depth.copy()
def write_raster(path, data, valid, profile, description, tags):
arr = np.where(valid, data, NODATA).astype(np.float32)
with rasterio.open(path, 'w', **profile) as dst:
dst.write(arr, 1)
dst.set_band_description(1, description)
dst.update_tags(**tags)
v = arr[valid].astype(np.float64)
dst.update_tags(1, STATISTICS_MINIMUM=str(v.min()),
STATISTICS_MAXIMUM=str(v.max()),
STATISTICS_MEAN=str(v.mean()),
STATISTICS_STDDEV=str(v.std()),
STATISTICS_VALID_PERCENT=str(100 * valid.mean()))
def save_results(prefix, depth, valid, profile, meta):
vals = depth[valid]
if not np.isfinite(vals).all():
raise ValueError('Non-finite predictions found')
lo, hi = np.percentile(vals, [1, 99])
if hi <= lo:
raise ValueError('Depth map is effectively constant; do not treat it as a valid result')
# Display only: linear sign reversal, NOT inverse-depth or metric height.
near = np.clip((hi - depth) / (hi - lo), 0, 1).astype(np.float32)
tags = dict(model='DA3MONO-LARGE', units='relative_arbitrary_units',
inference='whole_image', process_res=str(meta['process_res']),
native_prediction_shape=str(meta['native_prediction_shape']))
write_raster(str(prefix) + '_depth.tif', depth, valid, profile,
'Relative depth (larger = farther); NOT height in meters', tags)
write_raster(str(prefix) + '_near01.tif', near, valid, profile,
'Display only: near=1 far=0; clipped P1-P99; NOT height',
dict(tags, units='normalized_display_only',
display_formula='clip((P99-depth)/(P99-P1),0,1)',
display_p1=str(lo), display_p99=str(hi)))
preview = np.zeros(depth.shape, dtype=np.uint8)
preview[valid] = np.rint(near[valid] * 255).astype(np.uint8)
Image.fromarray(preview).save(str(prefix) + '_near.png')
meta.update(depth_statistics={
'min': float(vals.min()), 'max': float(vals.max()),
'p1': float(lo), 'median': float(np.median(vals)), 'p99': float(hi)},
valid_pixels=int(valid.sum()), nodata=NODATA,
depth_units='arbitrary; not meters',
preview='near=white; P1/P99 clipping; invalid=black',
output_resampling='bilinear onto original pixel grid; no new detail added')
Path(str(prefix) + '_metadata.json').write_text(
json.dumps(meta, ensure_ascii=False, indent=2), encoding='utf-8')
print(json.dumps(meta['depth_statistics'], indent=2), flush=True)
def main():
p = argparse.ArgumentParser(description=__doc__)
p.add_argument('--input', required=True)
p.add_argument('--output-prefix', required=True)
p.add_argument('--model', default=str(Path(__file__).resolve().parent / 'models/DA3MONO-LARGE'))
p.add_argument('--process-res', type=int, default=1008,
help='Inference longest side, divisible by 14; try 1008 then 2016')
p.add_argument('--bands', nargs=3, type=int, default=[1, 2, 3], metavar=('R', 'G', 'B'))
p.add_argument('--device', default='cuda:0')
args = p.parse_args()
if args.process_res < 28 or args.process_res % 14:
p.error('--process-res must be >=28 and divisible by 14')
inp = Path(args.input).expanduser().resolve()
prefix = Path(args.output_prefix).expanduser().resolve()
if not Path(args.model).is_dir():
p.error('Local model directory missing; run download_model.py first')
for suffix in ['_depth.tif', '_near01.tif', '_near.png', '_metadata.json']:
if Path(str(prefix) + suffix).exists():
p.error(f'Output exists: {prefix}{suffix}; use a new output prefix')
prefix.parent.mkdir(parents=True, exist_ok=True)
rgb, valid, profile = read_rgb(inp, args.bands)
height, width = valid.shape
if min(height, width) / max(height, width) * args.process_res < 28:
p.error('Image aspect ratio is too extreme at this process resolution')
print(f'Input: {width}x{height}, RGB bands {args.bands}; process-res={args.process_res}', flush=True)
import torch
from depth_anything_3.api import DepthAnything3
if not args.device.startswith('cuda') or not torch.cuda.is_available():
p.error('This workflow requires a CUDA GPU')
torch.cuda.set_device(torch.device(args.device))
print(f'GPU: {torch.cuda.get_device_name()}', flush=True)
model = DepthAnything3.from_pretrained(args.model).to(args.device).eval()
if model.model_name != 'da3mono-large':
raise ValueError(f'This workflow expects da3mono-large, got {model.model_name}')
start = time.time()
with torch.inference_mode():
prediction = model.inference(
[rgb], process_res=args.process_res,
process_res_method='upper_bound_resize',
infer_gs=False, export_dir=None,
)
native = np.asarray(prediction.depth[0], dtype=np.float32)
print('Native depth shape:', native.shape, flush=True)
depth = restore_depth(native, height, width)
meta = dict(input=str(inp), model=str(Path(args.model).resolve()),
model_name=model.model_name, process_res=args.process_res,
native_prediction_shape=list(native.shape), output_shape=[height, width],
crs=str(profile['crs']), transform=list(profile['transform']),
rgb_bands=args.bands, torch_version=torch.__version__,
cuda_runtime=torch.version.cuda, gpu=torch.cuda.get_device_name(),
inference_seconds=time.time() - start)
save_results(prefix, depth, valid, profile, meta)
print('Done:', prefix, flush=True)
if __name__ == '__main__':
main()
da3_tile_io.py
Full contents of da3_tile_io.py
"""Raster IO for experimental tiled relative depth mosaics."""
from pathlib import Path
import argparse
import json
import time
import numpy as np
import rasterio
from PIL import Image
NODATA = -9999.0
def read_rgb(path, bands):
with rasterio.open(path) as src:
if min(bands) < 1 or max(bands) > src.count:
raise ValueError(f'Invalid RGB bands {bands}; source has {src.count} bands')
rgb = src.read(bands)
if rgb.dtype != np.uint8:
raise ValueError('Input must be RGB uint8. Convert 16-bit/float data explicitly first.')
if src.crs is None:
raise ValueError('Input has no CRS; supply a georeferenced RGB GeoTIFF')
valid = np.all(src.read_masks(bands) > 0, axis=0)
if not valid.any():
raise ValueError('Input has no valid RGB pixels')
profile = dict(driver='GTiff', width=src.width, height=src.height,
count=1, dtype='float32', crs=src.crs, transform=src.transform,
nodata=NODATA, compress='deflate', predictor=3, tiled=True,
blockxsize=256, blockysize=256, BIGTIFF='IF_SAFER')
rgb = np.moveaxis(rgb, 0, -1).copy()
# Fill masked pixels only for inference; their outputs remain NoData.
if not valid.all():
rgb[~valid] = np.median(rgb[valid], axis=0).astype(np.uint8)
print('Warning: masked areas can affect nearby predictions.', flush=True)
return rgb, valid, profile
def restore_depth(depth, height, width):
depth = np.asarray(depth, dtype=np.float32)
if depth.ndim != 2 or not np.isfinite(depth).all():
raise ValueError('Model depth must be a finite 2-D array')
if depth.shape != (height, width):
depth = np.asarray(Image.fromarray(depth).resize(
(width, height), resample=Image.Resampling.BILINEAR), dtype=np.float32)
return depth.copy()
def write_raster(path, data, valid, profile, description, tags):
arr = np.where(valid, data, NODATA).astype(np.float32)
with rasterio.open(path, 'w', **profile) as dst:
dst.write(arr, 1)
dst.set_band_description(1, description)
dst.update_tags(**tags)
v = arr[valid].astype(np.float64)
dst.update_tags(1, STATISTICS_MINIMUM=str(v.min()),
STATISTICS_MAXIMUM=str(v.max()),
STATISTICS_MEAN=str(v.mean()),
STATISTICS_STDDEV=str(v.std()),
STATISTICS_VALID_PERCENT=str(100 * valid.mean()))
def save_results(prefix, depth, valid, profile, meta):
vals = depth[valid]
if not np.isfinite(vals).all():
raise ValueError('Non-finite predictions found')
lo, hi = np.percentile(vals, [1, 99])
if hi <= lo:
raise ValueError('Depth map is effectively constant; do not treat it as a valid result')
# Display only: linear sign reversal, NOT inverse-depth or metric height.
near = np.clip((hi - depth) / (hi - lo), 0, 1).astype(np.float32)
tags = dict(model='DA3MONO-LARGE', units='relative_arbitrary_units',
inference=meta['inference'], alignment=meta['alignment'], tile_size=str(meta['tile_size']), overlap=str(meta['overlap']), process_res=str(meta['process_res']),
native_prediction_shape=str(meta['native_prediction_shape']))
write_raster(str(prefix) + '_depth.tif', depth, valid, profile,
'Relative depth (larger = farther); NOT height in meters', tags)
write_raster(str(prefix) + '_near01.tif', near, valid, profile,
'Display only: near=1 far=0; clipped P1-P99; NOT height',
dict(tags, units='normalized_display_only',
display_formula='clip((P99-depth)/(P99-P1),0,1)',
display_p1=str(lo), display_p99=str(hi)))
preview = np.zeros(depth.shape, dtype=np.uint8)
preview[valid] = np.rint(near[valid] * 255).astype(np.uint8)
Image.fromarray(preview).save(str(prefix) + '_near.png')
meta.update(depth_statistics={
'min': float(vals.min()), 'max': float(vals.max()),
'p1': float(lo), 'median': float(np.median(vals)), 'p99': float(hi)},
valid_pixels=int(valid.sum()), nodata=NODATA,
depth_units='arbitrary; not meters',
preview='near=white; P1/P99 clipping; invalid=black',
output_resampling='bilinear onto original pixel grid; no new detail added')
Path(str(prefix) + '_metadata.json').write_text(
json.dumps(meta, ensure_ascii=False, indent=2), encoding='utf-8')
print(json.dumps(meta['depth_statistics'], indent=2), flush=True)
infer_geotiff_tiled.py
Full contents of infer_geotiff_tiled.py
"""DA3MONO independent overlapping tile inference; experimental relative-depth mosaic."""
import argparse
from pathlib import Path
import time
import numpy as np
from da3_tile_io import read_rgb, restore_depth, save_results
def starts(length, tile, overlap):
if length <= tile:
return [0]
pos = list(range(0, length - tile + 1, tile - overlap))
if pos[-1] != length - tile:
pos.append(length - tile)
return pos
def fit_overlap(x, y, mode):
"""Robust positive affine fit; reject poorly constrained scale estimates."""
x, y = np.asarray(x, dtype=np.float64), np.asarray(y, dtype=np.float64)
if not len(x):
return 1., 0., 'unanchored'
step = max(1, len(x) // 20000)
x, y = x[::step], y[::step]
shift = float(np.median(y - x))
if mode == 'shift' or len(x) < 128:
return 1., shift, 'shift'
keep = np.ones(len(x), dtype=bool)
a, b = 1., shift
for _ in range(4):
xx, yy = x[keep], y[keep]
if len(xx) < 128 or np.std(xx) < 1e-6 * max(1., abs(np.mean(xx))):
return 1., shift, 'fallback_flat_overlap'
dx, dy = xx - xx.mean(), yy - yy.mean()
corr = np.dot(dx, dy) / max(np.linalg.norm(dx) * np.linalg.norm(dy), 1e-20)
a = float(np.dot(dx, dy) / np.dot(dx, dx))
if corr < .5 or not .25 <= a <= 4.:
return 1., shift, 'fallback_unstable_scale'
b = float(np.median(yy - a * xx))
residual = y - (a * x + b)
med = np.median(residual)
mad = np.median(np.abs(residual - med))
keep = np.abs(residual - med) <= max(3 * 1.4826 * mad, 1e-6)
return a, b, 'affine'
def blend_weight(height, width):
# Positive edge weights retain image boundary pixels.
wy = np.maximum(np.sin(np.pi * (np.arange(height) + .5) / height)**2, .02)
wx = np.maximum(np.sin(np.pi * (np.arange(width) + .5) / width)**2, .02)
return wy[:, None] * wx[None, :]
def mosaic(rgb, valid, predict, tile, overlap, alignment):
h, w = valid.shape
total = np.zeros((h, w), dtype=np.float64)
weights = np.zeros((h, w), dtype=np.float64)
records = []
ys, xs = starts(h, tile, overlap), starts(w, tile, overlap)
for y in ys:
for x in xs:
sl = np.s_[y:min(y+tile, h), x:min(x+tile, w)]
mask = valid[sl]
if not mask.any():
continue
patch = rgb[sl]
ph, pw = mask.shape
# Small images padded; normal edge tiles remain exactly tile x tile.
padded = np.pad(patch, ((0, tile-ph), (0, tile-pw), (0, 0)), mode='edge')
native = predict(padded)
depth = restore_depth(native, tile, tile)[:ph, :pw]
old_w, old_sum = weights[sl], total[sl]
common = mask & (old_w > 0)
a, b, status = 1., 0., 'none'
if alignment != 'none':
a, b, status = fit_overlap(depth[common], old_sum[common]/old_w[common], alignment)
corrected = a * depth + b
weight = blend_weight(ph, pw) * mask
old_sum += corrected * weight
old_w += weight
records.append(dict(x=x, y=y, width=pw, height=ph,
native_shape=list(np.shape(native)), scale=a, shift=b,
alignment_status=status, overlap_pixels=int(common.sum())))
print(f'Tile {len(records)}/{len(ys)*len(xs)} ({x},{y}) {status}: a={a:.4g}, b={b:.4g}', flush=True)
if np.any(valid & (weights <= 0)):
raise RuntimeError('Valid pixels not covered')
output = np.zeros((h, w), dtype=np.float32)
output[valid] = (total[valid] / weights[valid]).astype(np.float32)
return output, records
def main():
p = argparse.ArgumentParser(description=__doc__)
p.add_argument('--input', required=True)
p.add_argument('--output-prefix', required=True)
p.add_argument('--model', default=str(Path(__file__).resolve().parent/'models/DA3MONO-LARGE'))
p.add_argument('--tile-size', type=int, default=300)
p.add_argument('--overlap', type=int, default=100, help='Overlap pixels, NOT stride')
p.add_argument('--process-res', type=int, default=518, help='Model input longest side, multiple of 14')
p.add_argument('--alignment', choices=['affine', 'shift', 'none'], default='affine')
p.add_argument('--bands', nargs=3, type=int, default=[1,2,3])
p.add_argument('--device', default='cuda:0')
args = p.parse_args()
if args.tile_size < 28 or not 0 <= args.overlap < args.tile_size:
p.error('Require tile-size >=28 and 0 <= overlap < tile-size')
if args.alignment != 'none' and args.overlap == 0:
p.error('Alignment needs overlap >0; use --alignment none for nonoverlapping tiles')
if args.process_res < 28 or args.process_res % 14:
p.error('process-res must be >=28 and divisible by 14')
if not Path(args.model).is_dir():
p.error('Model directory missing; specify existing --model path')
prefix = Path(args.output_prefix).expanduser().resolve()
for suffix in ['_depth.tif','_near01.tif','_near.png','_metadata.json']:
if Path(str(prefix)+suffix).exists():
p.error(f'Output exists: {prefix}{suffix}; use a new prefix')
rgb, valid, profile = read_rgb(args.input, args.bands)
import torch
from depth_anything_3.api import DepthAnything3
if not args.device.startswith('cuda') or not torch.cuda.is_available():
p.error('CUDA GPU required')
torch.cuda.set_device(torch.device(args.device))
model = DepthAnything3.from_pretrained(args.model).to(args.device).eval()
if model.model_name != 'da3mono-large':
p.error('Expected DA3MONO-LARGE')
def predict(patch):
with torch.inference_mode():
result = model.inference([patch], process_res=args.process_res,
process_res_method='upper_bound_resize',
infer_gs=False, export_dir=None)
return np.asarray(result.depth[0], dtype=np.float32)
start = time.time()
depth, records = mosaic(rgb, valid, predict, args.tile_size, args.overlap, args.alignment)
meta = dict(input=str(Path(args.input).resolve()), model=str(Path(args.model).resolve()),
model_name=model.model_name, process_res=args.process_res,
native_prediction_shape=records[0]['native_shape'], output_shape=list(valid.shape),
crs=str(profile['crs']), transform=list(profile['transform']), rgb_bands=args.bands,
inference='independent_overlapping_tiles', tile_size=args.tile_size,
overlap=args.overlap, stride=args.tile_size-args.overlap, alignment=args.alignment,
alignment_warning='Experimental relative mosaic; scan-order drift and disconnected scales possible; not metric depth or height',
tiles=records, torch_version=torch.__version__, cuda_runtime=torch.version.cuda,
gpu=torch.cuda.get_device_name(), inference_seconds=time.time()-start)
prefix.parent.mkdir(parents=True, exist_ok=True)
save_results(prefix, depth, valid, profile, meta)
print('Done:', prefix, flush=True)
if __name__ == '__main__':
main()
Top comments (0)