Skip to content

ds_msp.models

The seven camera models, each implementing the CameraModel contract. Pick by FOV and required accuracy — see Choosing a camera model.

Double Sphere

ds_msp.models.DoubleSphereModel

Double Sphere camera (Usenko et al. 2018). Satisfies CameraModel.

Source code in ds_msp/models/double_sphere.py 
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
class DoubleSphereModel:
    """Double Sphere camera (Usenko et al. 2018). Satisfies ``CameraModel``."""

    name: ClassVar[str] = "ds"
    param_names: ClassVar[Tuple[str, ...]] = ("fx", "fy", "cx", "cy", "xi", "alpha")

    def __init__(self, fx: float, fy: float, cx: float, cy: float,
                 xi: float, alpha: float) -> None:
        self.fx = float(fx)
        self.fy = float(fy)
        self.cx = float(cx)
        self.cy = float(cy)
        self.xi = float(xi)
        self.alpha = float(alpha)

    @classmethod
    def sample(cls) -> "DoubleSphereModel":
        """Realistic instance for contract testing (the bundled calibration)."""
        return cls(711.57, 711.24, 949.18, 518.81, 0.183, 0.809)

    # -- parameter access -------------------------------------------------
    @property
    def params(self) -> np.ndarray:
        """Flat parameter vector ``[fx, fy, cx, cy, xi, alpha]``, see `param_names`."""
        return np.array([self.fx, self.fy, self.cx, self.cy, self.xi, self.alpha],
                        dtype=np.float64)

    @property
    def K(self) -> np.ndarray:
        """3x3 pinhole intrinsic matrix built from ``fx, fy, cx, cy``."""
        return np.array([[self.fx, 0.0, self.cx],
                         [0.0, self.fy, self.cy],
                         [0.0, 0.0, 1.0]], dtype=np.float64)

    @property
    def distortion(self) -> np.ndarray:
        """Distortion tail ``[xi, alpha]`` (sphere offset, perspective blend)."""
        return np.array([self.xi, self.alpha], dtype=np.float64)

    # -- core math --------------------------------------------------------
    def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Project camera-frame points to pixels via the two-sphere composition.

        The Double Sphere model is what makes this model's family distinct: a
        point is first re-centered onto a unit sphere shifted by ``xi`` along
        -Z (``z1 = z + xi*d1``), then perspective-divided from a second point
        blended between that shifted sphere and the pinhole plane by
        ``alpha``. Composing two spheres (rather than one, as in `UCMModel`)
        is what lets this closed-form model reach fields of view beyond 180°
        (Usenko et al. 2018, Sec. 3).

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates ``(u, v)``, origin top-left.
        valid : ndarray, shape (...,)
            Projectability mask. ``True`` iff the point lies in the tilted
            half-space ``z > -w2 * d1`` (*not* the naive ``z > 0``) and the
            projection denominator is bounded away from zero — see
            ``ds_msp.models.ds_math.ds_project`` and
            [Projection validity and FOV](https://github.com/Munna-Manoj/DS-MSP/blob/main/docs/explain/projection_validity_and_fov.md).
            Invalid entries are masked, never NaN.

        References
        ----------
        Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model."
        3DV 2018 (projection; see Eq. 43-45 for the validity half-space).

        Examples
        --------
        A point *behind* the pinhole plane (``z < 0``) can still be valid —
        this is the wide-FOV behavior the two-sphere composition exists for:

        >>> import numpy as np
        >>> from ds_msp.models import DoubleSphereModel
        >>> m = DoubleSphereModel.sample()
        >>> uv, valid = m.project(np.array([[1.0, 0.0, -0.3]]))
        >>> bool(valid[0])
        True
        """
        u, v, valid = ds_project(np.asarray(P, dtype=np.float64),
                                 self.fx, self.fy, self.cx, self.cy, self.xi, self.alpha)
        return np.stack([u, v], axis=-1), valid

    def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Unproject pixels to unit bearing rays (closed form).

        Inverts the two-sphere composition analytically (no iteration): first
        recovers the point on the shifted sphere from the normalized pixel,
        then un-shifts by ``xi``. See ``ds_msp.models.ds_math.ds_unproject``
        for the exact algebra.

        Parameters
        ----------
        uv : ndarray, shape (..., 2)
            Pixel coordinates ``(u, v)``, origin top-left.

        Returns
        -------
        rays : ndarray, shape (..., 3)
            Unit-norm camera-frame bearing vectors, +Z forward.
        valid : ndarray, shape (...,)
            ``True`` iff the pixel lies within the sphere's reachable disc
            (radius set by ``alpha``) and every intermediate square root is
            real and the perspective denominator is bounded away from zero.
            Invalid rays are zeroed, never NaN.

        References
        ----------
        Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model." 3DV 2018
        (closed-form unprojection).

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import DoubleSphereModel
        >>> m = DoubleSphereModel.sample()
        >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
        >>> np.round(rays, 4)
        array([[0., 0., 1.]])
        """
        return ds_unproject(np.asarray(uv, dtype=np.float64),
                            self.fx, self.fy, self.cx, self.cy, self.xi, self.alpha)

    def project_jacobian(
        self, P: np.ndarray
    ) -> Tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray]:
        """Project with analytic derivatives (no autodiff, no finite differences).

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates, identical to `project`.
        J_point : ndarray, shape (..., 2, 3)
            ``d(u, v) / d(x, y, z)``.
        J_param : ndarray, shape (..., 2, 6)
            ``d(u, v) / d(fx, fy, cx, cy, xi, alpha)``, columns in
            `param_names` order.
        valid : ndarray, shape (...,)
            Projectability mask, identical condition to `project`.

        References
        ----------
        Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model." 3DV 2018
        (closed-form Jacobian derived from the projection above; verified here by
        finite-difference check, relative error <= 1e-6, see ``pytest -m jac``).

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import DoubleSphereModel
        >>> m = DoubleSphereModel.sample()
        >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
        >>> J_point.shape, J_param.shape
        ((1, 2, 3), (1, 2, 6))
        """
        u, v, J_point, J_param, valid = ds_project_jacobian(
            np.asarray(P, dtype=np.float64),
            self.fx, self.fy, self.cx, self.cy, self.xi, self.alpha)
        return np.stack([u, v], axis=-1), J_point, J_param, valid

    # -- construction / bounds -------------------------------------------
    @classmethod
    def from_params(cls, p: np.ndarray) -> "DoubleSphereModel":
        """Build from a flat ``[fx, fy, cx, cy, xi, alpha]`` vector."""
        return cls(*np.asarray(p, dtype=np.float64).ravel())

    @classmethod
    def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
        """Optimizer bounds: ``xi in [-1, 1]``, ``alpha in (0, 1)``, focal/center wide-open."""
        lb = np.array([1.0, 1.0, -1e5, -1e5, -1.0, 1e-6], dtype=np.float64)
        ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0, 1.0 - 1e-6], dtype=np.float64)
        return lb, ub

    # -- conversion hook --------------------------------------------------
    def initialize_from_correspondences(
        self, K_seed: np.ndarray, rays: np.ndarray, pixels: np.ndarray
    ) -> None:
        """Seed ``fx,fy,cx,cy`` from `K_seed`; seed ``xi=0`` and solve ``alpha`` linearly
        (reduces to the UCM closed form since ``xi=0`` collapses both spheres into one)."""
        # Inherit pinhole intrinsics from the source.
        self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
        self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
        # Seed distortion: xi = 0 reduces DS to UCM; for unit rays (d1 = 1),
        #   mx = x / (alpha*(1 - z) + z)  =>  alpha = (x - mx*z) / (mx*(1 - z)).
        # Solve linearly over both axes.
        rays = np.asarray(rays, dtype=np.float64)
        x, y, z = rays[:, 0], rays[:, 1], rays[:, 2]
        mx = (pixels[:, 0] - self.cx) / self.fx
        my = (pixels[:, 1] - self.cy) / self.fy
        A = np.concatenate([mx * (1.0 - z), my * (1.0 - z)])
        b = np.concatenate([x - mx * z, y - my * z])
        denom = float(A @ A)
        self.xi = 0.0
        self.alpha = float(np.clip((A @ b) / denom, 1e-6, 1.0 - 1e-6)) if denom > 1e-12 else 0.5

    # -- serialization ----------------------------------------------------
    def to_dict(self) -> dict:
        """Serialize to ``{"model": "ds", "fx": ..., ..., "alpha": ...}``."""
        d = {"model": self.name}
        d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
        return d

    @classmethod
    def from_dict(cls, d: dict) -> "DoubleSphereModel":
        """Reconstruct from :meth:`to_dict` output."""
        return cls(**{k: d[k] for k in cls.param_names})

    def __repr__(self) -> str:
        return ("DoubleSphereModel(fx={:.3f}, fy={:.3f}, cx={:.3f}, cy={:.3f}, "
                "xi={:.4f}, alpha={:.4f})").format(
                    self.fx, self.fy, self.cx, self.cy, self.xi, self.alpha)

K property 

K: ndarray

3x3 pinhole intrinsic matrix built from fx, fy, cx, cy.

distortion property 

distortion: ndarray

Distortion tail [xi, alpha] (sphere offset, perspective blend).

params property 

params: ndarray

Flat parameter vector [fx, fy, cx, cy, xi, alpha], see param_names.

from_dict classmethod 

from_dict(d: dict) -> 'DoubleSphereModel'

Reconstruct from :meth:to_dict output.

Source code in ds_msp/models/double_sphere.py 
230
231
232
233
@classmethod
def from_dict(cls, d: dict) -> "DoubleSphereModel":
    """Reconstruct from :meth:`to_dict` output."""
    return cls(**{k: d[k] for k in cls.param_names})

from_params classmethod 

from_params(p: ndarray) -> 'DoubleSphereModel'

Build from a flat [fx, fy, cx, cy, xi, alpha] vector.

Source code in ds_msp/models/double_sphere.py 
189
190
191
192
@classmethod
def from_params(cls, p: np.ndarray) -> "DoubleSphereModel":
    """Build from a flat ``[fx, fy, cx, cy, xi, alpha]`` vector."""
    return cls(*np.asarray(p, dtype=np.float64).ravel())

initialize_from_correspondences 

initialize_from_correspondences(K_seed: ndarray, rays: ndarray, pixels: ndarray) -> None

Seed fx,fy,cx,cy from K_seed; seed xi=0 and solve alpha linearly (reduces to the UCM closed form since xi=0 collapses both spheres into one).

Source code in ds_msp/models/double_sphere.py 
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
def initialize_from_correspondences(
    self, K_seed: np.ndarray, rays: np.ndarray, pixels: np.ndarray
) -> None:
    """Seed ``fx,fy,cx,cy`` from `K_seed`; seed ``xi=0`` and solve ``alpha`` linearly
    (reduces to the UCM closed form since ``xi=0`` collapses both spheres into one)."""
    # Inherit pinhole intrinsics from the source.
    self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
    self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
    # Seed distortion: xi = 0 reduces DS to UCM; for unit rays (d1 = 1),
    #   mx = x / (alpha*(1 - z) + z)  =>  alpha = (x - mx*z) / (mx*(1 - z)).
    # Solve linearly over both axes.
    rays = np.asarray(rays, dtype=np.float64)
    x, y, z = rays[:, 0], rays[:, 1], rays[:, 2]
    mx = (pixels[:, 0] - self.cx) / self.fx
    my = (pixels[:, 1] - self.cy) / self.fy
    A = np.concatenate([mx * (1.0 - z), my * (1.0 - z)])
    b = np.concatenate([x - mx * z, y - my * z])
    denom = float(A @ A)
    self.xi = 0.0
    self.alpha = float(np.clip((A @ b) / denom, 1e-6, 1.0 - 1e-6)) if denom > 1e-12 else 0.5

param_bounds classmethod 

param_bounds() -> Tuple[np.ndarray, np.ndarray]

Optimizer bounds: xi in [-1, 1], alpha in (0, 1), focal/center wide-open.

Source code in ds_msp/models/double_sphere.py 
194
195
196
197
198
199
@classmethod
def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
    """Optimizer bounds: ``xi in [-1, 1]``, ``alpha in (0, 1)``, focal/center wide-open."""
    lb = np.array([1.0, 1.0, -1e5, -1e5, -1.0, 1e-6], dtype=np.float64)
    ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0, 1.0 - 1e-6], dtype=np.float64)
    return lb, ub

project 

project(P: ndarray) -> Tuple[np.ndarray, np.ndarray]

Project camera-frame points to pixels via the two-sphere composition.

The Double Sphere model is what makes this model's family distinct: a point is first re-centered onto a unit sphere shifted by xi along -Z (z1 = z + xi*d1), then perspective-divided from a second point blended between that shifted sphere and the pinhole plane by alpha. Composing two spheres (rather than one, as in UCMModel) is what lets this closed-form model reach fields of view beyond 180° (Usenko et al. 2018, Sec. 3).

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates (u, v), origin top-left.
valid (ndarray, shape(...))

Projectability mask. True iff the point lies in the tilted half-space z > -w2 * d1 (not the naive z > 0) and the projection denominator is bounded away from zero — see ds_msp.models.ds_math.ds_project and Projection validity and FOV. Invalid entries are masked, never NaN.
References

Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model." 3DV 2018 (projection; see Eq. 43-45 for the validity half-space).

Examples:

A point behind the pinhole plane (z < 0) can still be valid — this is the wide-FOV behavior the two-sphere composition exists for:

>>> import numpy as np
>>> from ds_msp.models import DoubleSphereModel
>>> m = DoubleSphereModel.sample()
>>> uv, valid = m.project(np.array([[1.0, 0.0, -0.3]]))
>>> bool(valid[0])
True
Source code in ds_msp/models/double_sphere.py 
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Project camera-frame points to pixels via the two-sphere composition.

    The Double Sphere model is what makes this model's family distinct: a
    point is first re-centered onto a unit sphere shifted by ``xi`` along
    -Z (``z1 = z + xi*d1``), then perspective-divided from a second point
    blended between that shifted sphere and the pinhole plane by
    ``alpha``. Composing two spheres (rather than one, as in `UCMModel`)
    is what lets this closed-form model reach fields of view beyond 180°
    (Usenko et al. 2018, Sec. 3).

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates ``(u, v)``, origin top-left.
    valid : ndarray, shape (...,)
        Projectability mask. ``True`` iff the point lies in the tilted
        half-space ``z > -w2 * d1`` (*not* the naive ``z > 0``) and the
        projection denominator is bounded away from zero — see
        ``ds_msp.models.ds_math.ds_project`` and
        [Projection validity and FOV](https://github.com/Munna-Manoj/DS-MSP/blob/main/docs/explain/projection_validity_and_fov.md).
        Invalid entries are masked, never NaN.

    References
    ----------
    Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model."
    3DV 2018 (projection; see Eq. 43-45 for the validity half-space).

    Examples
    --------
    A point *behind* the pinhole plane (``z < 0``) can still be valid —
    this is the wide-FOV behavior the two-sphere composition exists for:

    >>> import numpy as np
    >>> from ds_msp.models import DoubleSphereModel
    >>> m = DoubleSphereModel.sample()
    >>> uv, valid = m.project(np.array([[1.0, 0.0, -0.3]]))
    >>> bool(valid[0])
    True
    """
    u, v, valid = ds_project(np.asarray(P, dtype=np.float64),
                             self.fx, self.fy, self.cx, self.cy, self.xi, self.alpha)
    return np.stack([u, v], axis=-1), valid

project_jacobian 

project_jacobian(P: ndarray) -> Tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray]

Project with analytic derivatives (no autodiff, no finite differences).

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates, identical to project.
J_point (ndarray, shape(..., 2, 3))

d(u, v) / d(x, y, z).
J_param (ndarray, shape(..., 2, 6))

d(u, v) / d(fx, fy, cx, cy, xi, alpha), columns in param_names order.
valid (ndarray, shape(...))

Projectability mask, identical condition to project.
References

Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model." 3DV 2018 (closed-form Jacobian derived from the projection above; verified here by finite-difference check, relative error <= 1e-6, see pytest -m jac).

Examples:

>>> import numpy as np
>>> from ds_msp.models import DoubleSphereModel
>>> m = DoubleSphereModel.sample()
>>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
>>> J_point.shape, J_param.shape
((1, 2, 3), (1, 2, 6))
Source code in ds_msp/models/double_sphere.py 
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
def project_jacobian(
    self, P: np.ndarray
) -> Tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray]:
    """Project with analytic derivatives (no autodiff, no finite differences).

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates, identical to `project`.
    J_point : ndarray, shape (..., 2, 3)
        ``d(u, v) / d(x, y, z)``.
    J_param : ndarray, shape (..., 2, 6)
        ``d(u, v) / d(fx, fy, cx, cy, xi, alpha)``, columns in
        `param_names` order.
    valid : ndarray, shape (...,)
        Projectability mask, identical condition to `project`.

    References
    ----------
    Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model." 3DV 2018
    (closed-form Jacobian derived from the projection above; verified here by
    finite-difference check, relative error <= 1e-6, see ``pytest -m jac``).

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import DoubleSphereModel
    >>> m = DoubleSphereModel.sample()
    >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
    >>> J_point.shape, J_param.shape
    ((1, 2, 3), (1, 2, 6))
    """
    u, v, J_point, J_param, valid = ds_project_jacobian(
        np.asarray(P, dtype=np.float64),
        self.fx, self.fy, self.cx, self.cy, self.xi, self.alpha)
    return np.stack([u, v], axis=-1), J_point, J_param, valid

sample classmethod 

sample() -> 'DoubleSphereModel'

Realistic instance for contract testing (the bundled calibration).

Source code in ds_msp/models/double_sphere.py 
32
33
34
35
@classmethod
def sample(cls) -> "DoubleSphereModel":
    """Realistic instance for contract testing (the bundled calibration)."""
    return cls(711.57, 711.24, 949.18, 518.81, 0.183, 0.809)

to_dict 

to_dict() -> dict

Serialize to {"model": "ds", "fx": ..., ..., "alpha": ...}.

Source code in ds_msp/models/double_sphere.py 
224
225
226
227
228
def to_dict(self) -> dict:
    """Serialize to ``{"model": "ds", "fx": ..., ..., "alpha": ...}``."""
    d = {"model": self.name}
    d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
    return d

unproject 

unproject(uv: ndarray) -> Tuple[np.ndarray, np.ndarray]

Unproject pixels to unit bearing rays (closed form).

Inverts the two-sphere composition analytically (no iteration): first recovers the point on the shifted sphere from the normalized pixel, then un-shifts by xi. See ds_msp.models.ds_math.ds_unproject for the exact algebra.

Parameters:

Name Type Description Default
uv (ndarray, shape(..., 2))

Pixel coordinates (u, v), origin top-left.

 required

Returns:

Name Type Description
rays (ndarray, shape(..., 3))

Unit-norm camera-frame bearing vectors, +Z forward.
valid (ndarray, shape(...))

True iff the pixel lies within the sphere's reachable disc (radius set by alpha) and every intermediate square root is real and the perspective denominator is bounded away from zero. Invalid rays are zeroed, never NaN.
References

Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model." 3DV 2018 (closed-form unprojection).

Examples:

>>> import numpy as np
>>> from ds_msp.models import DoubleSphereModel
>>> m = DoubleSphereModel.sample()
>>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
>>> np.round(rays, 4)
array([[0., 0., 1.]])
Source code in ds_msp/models/double_sphere.py 
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Unproject pixels to unit bearing rays (closed form).

    Inverts the two-sphere composition analytically (no iteration): first
    recovers the point on the shifted sphere from the normalized pixel,
    then un-shifts by ``xi``. See ``ds_msp.models.ds_math.ds_unproject``
    for the exact algebra.

    Parameters
    ----------
    uv : ndarray, shape (..., 2)
        Pixel coordinates ``(u, v)``, origin top-left.

    Returns
    -------
    rays : ndarray, shape (..., 3)
        Unit-norm camera-frame bearing vectors, +Z forward.
    valid : ndarray, shape (...,)
        ``True`` iff the pixel lies within the sphere's reachable disc
        (radius set by ``alpha``) and every intermediate square root is
        real and the perspective denominator is bounded away from zero.
        Invalid rays are zeroed, never NaN.

    References
    ----------
    Usenko, V., Demmel, N., Cremers, D. "The Double Sphere Camera Model." 3DV 2018
    (closed-form unprojection).

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import DoubleSphereModel
    >>> m = DoubleSphereModel.sample()
    >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
    >>> np.round(rays, 4)
    array([[0., 0., 1.]])
    """
    return ds_unproject(np.asarray(uv, dtype=np.float64),
                        self.fx, self.fy, self.cx, self.cy, self.xi, self.alpha)

Double Sphere⁺ (DS⁺)

ds_msp.models.DSPlusModel

DS+ (UCM core + division radial + 2-axis tilt). Satisfies CameraModel.

Source code in ds_msp/models/dsplus.py 
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
class DSPlusModel:
    """DS+ (UCM core + division radial + 2-axis tilt). Satisfies ``CameraModel``."""

    name: ClassVar[str] = "dsplus"
    param_names: ClassVar[Tuple[str, ...]] = (
        "fx", "fy", "cx", "cy", "alpha", "lambda1", "lambda2", "tau_x", "tau_y")

    def __init__(self, fx, fy, cx, cy, alpha=0.5,
                 lambda1=0.0, lambda2=0.0, tau_x=0.0, tau_y=0.0) -> None:
        self.fx = float(fx)
        self.fy = float(fy)
        self.cx = float(cx)
        self.cy = float(cy)
        self.alpha = float(alpha)
        self.lambda1 = float(lambda1)
        self.lambda2 = float(lambda2)
        self.tau_x = float(tau_x)
        self.tau_y = float(tau_y)

    @classmethod
    def sample(cls) -> "DSPlusModel":
        """Realistic instance for contract testing (the bundled calibration)."""
        return cls(711.57, 711.24, 949.18, 518.81, 0.62, -0.10, 0.02, 0.001, -0.001)

    @property
    def params(self) -> np.ndarray:
        """Flat parameter vector in `param_names` order."""
        return np.array([self.fx, self.fy, self.cx, self.cy, self.alpha,
                         self.lambda1, self.lambda2, self.tau_x, self.tau_y],
                        dtype=np.float64)

    @property
    def K(self) -> np.ndarray:
        """3x3 pinhole intrinsic matrix built from ``fx, fy, cx, cy``."""
        return np.array([[self.fx, 0.0, self.cx], [0.0, self.fy, self.cy],
                         [0.0, 0.0, 1.0]], dtype=np.float64)

    @property
    def distortion(self) -> np.ndarray:
        """Distortion tail ``[alpha, lambda1, lambda2, tau_x, tau_y]``."""
        return np.array([self.alpha, self.lambda1, self.lambda2,
                         self.tau_x, self.tau_y], dtype=np.float64)

    def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Project camera-frame points through 4 closed-form-invertible stages.

        DS+ replaces Double Sphere's ``xi`` shift (the proven near-null
        direction for the target lens, see the module docstring) with a UCM
        sphere core followed by two extra closed-form stages: a Fitzgibbon
        division-model radial layer and a 2-axis Scheimpflug tilt homography,
        so the composed map stays fully closed-form invertible end to end::

            pixel = K . H_tau . D_lambda . S_alpha(bearing)

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates ``(u, v)``, origin top-left.
        valid : ndarray, shape (...,)
            ``True`` iff all 4 stages are individually valid: the UCM
            half-space/denominator condition, a positive division-model
            radial factor, and a nonzero tilt-homography denominator. See
            ``ds_msp.models.dsplus_math.dsplus_project``.

        References
        ----------
        UCM core: Geyer, C., Daniilidis, K. "A Unifying Theory for Central
        Panoramic Systems." ECCV 2000; Mei, C., Rives, P. ICRA 2007.
        Division radial layer: Fitzgibbon, A. "Simultaneous linear estimation
        of multiple view geometry and lens distortion." CVPR 2001.
        Tilt homography: cf. OpenCV ``CALIB_TILTED_MODEL``. Staged
        composition and Jacobian are this repo's own extension — see
        ``ds_msp/models/dsplus_math.py`` and
        [ADR-0005](https://github.com/Munna-Manoj/DS-MSP/blob/main/docs/process/architecture/decisions/ADR-0005-dsplus-eucmplus.md).

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import DSPlusModel
        >>> m = DSPlusModel.sample()
        >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
        >>> np.round(uv, 2)
        array([[949.18, 518.81]])
        """
        u, v, valid = dsplus_project(
            np.asarray(P, dtype=np.float64),
            self.fx, self.fy, self.cx, self.cy, self.alpha,
            self.lambda1, self.lambda2, self.tau_x, self.tau_y)
        return np.stack([u, v], axis=-1), valid

    def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Unproject pixels to unit bearing rays (closed form, no iteration).

        Inverts the 4-stage chain in reverse (``K^-1``, tilt inverse, then a
        quadratic/quartic radical to invert the division-model radial layer,
        then the UCM closed form). See
        ``ds_msp.models.dsplus_math.dsplus_unproject`` and
        ``ds_msp.models.dsplus_math._invert_division``.

        Parameters
        ----------
        uv : ndarray, shape (..., 2)
            Pixel coordinates ``(u, v)``, origin top-left.

        Returns
        -------
        rays : ndarray, shape (..., 3)
            Unit-norm camera-frame bearing vectors, +Z forward.
        valid : ndarray, shape (...,)
            ``True`` iff the tilt inverse, the division-radial root, and the
            UCM sphere inverse are all well-defined. Invalid rays are
            zeroed, never NaN.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import DSPlusModel
        >>> m = DSPlusModel.sample()
        >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
        >>> np.round(rays, 4)
        array([[0., 0., 1.]])
        """
        return dsplus_unproject(
            np.asarray(uv, dtype=np.float64),
            self.fx, self.fy, self.cx, self.cy, self.alpha,
            self.lambda1, self.lambda2, self.tau_x, self.tau_y)

    def project_jacobian(self, P):
        """Project with analytic derivatives via the chain rule through all 4 stages.

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates, identical to `project`.
        J_point : ndarray, shape (..., 2, 3)
            ``d(u, v) / d(x, y, z)``, computed as ``K @ J_H @ J_D @ J_S``.
        J_param : ndarray, shape (..., 2, 9)
            ``d(u, v) / d(fx, fy, cx, cy, alpha, lambda1, lambda2, tau_x, tau_y)``,
            columns in `param_names` order.
        valid : ndarray, shape (...,)
            Projectability mask, identical condition to `project`.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import DSPlusModel
        >>> m = DSPlusModel.sample()
        >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
        >>> J_point.shape, J_param.shape
        ((1, 2, 3), (1, 2, 9))
        """
        u, v, J_point, J_param, valid = dsplus_project_jacobian(
            np.asarray(P, dtype=np.float64),
            self.fx, self.fy, self.cx, self.cy, self.alpha,
            self.lambda1, self.lambda2, self.tau_x, self.tau_y)
        return np.stack([u, v], axis=-1), J_point, J_param, valid

    @classmethod
    def from_params(cls, p: np.ndarray) -> "DSPlusModel":
        """Build from a flat vector in `param_names` order."""
        return cls(*np.asarray(p, dtype=np.float64).ravel())

    @classmethod
    def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
        """Optimizer bounds: ``alpha in (0, 1)``, tilt/radial terms modestly bounded."""
        lb = np.array([1.0, 1.0, -1e5, -1e5, 1e-6, -2.0, -2.0, -0.5, -0.5],
                      dtype=np.float64)
        ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0 - 1e-6, 2.0, 2.0, 0.5, 0.5],
                      dtype=np.float64)
        return lb, ub

    def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
        """Seed ``fx,fy,cx,cy`` from `K_seed`; solve ``alpha`` linearly (same UCM
        closed form as `UCMModel`); zero the radial/tilt terms."""
        self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
        self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
        rays = np.asarray(rays, dtype=np.float64)
        x, y, z = rays[:, 0], rays[:, 1], rays[:, 2]
        mx = (pixels[:, 0] - self.cx) / self.fx
        my = (pixels[:, 1] - self.cy) / self.fy
        # Same UCM linear solve as ucm.py: alpha = (x - mx*z)/(mx*(1-z)).
        A = np.concatenate([mx * (1.0 - z), my * (1.0 - z)])
        b = np.concatenate([x - mx * z, y - my * z])
        denom = float(A @ A)
        self.alpha = float(np.clip((A @ b) / denom, 1e-6, 1.0 - 1e-6)) if denom > 1e-12 else 0.5
        self.lambda1 = 0.0
        self.lambda2 = 0.0
        self.tau_x = 0.0
        self.tau_y = 0.0

    def to_dict(self) -> dict:
        """Serialize to ``{"model": "dsplus", "fx": ..., ..., "tau_y": ...}``."""
        d = {"model": self.name}
        d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
        return d

    @classmethod
    def from_dict(cls, d: dict) -> "DSPlusModel":
        """Reconstruct from :meth:`to_dict` output."""
        return cls(**{k: d[k] for k in cls.param_names})

    def __repr__(self) -> str:
        return ("DSPlusModel(fx={:.3f}, fy={:.3f}, cx={:.3f}, cy={:.3f}, "
                "alpha={:.4f}, lambda1={:.5f}, lambda2={:.5f}, "
                "tau_x={:.5f}, tau_y={:.5f})").format(
                    self.fx, self.fy, self.cx, self.cy, self.alpha,
                    self.lambda1, self.lambda2, self.tau_x, self.tau_y)

K property 

K: ndarray

3x3 pinhole intrinsic matrix built from fx, fy, cx, cy.

distortion property 

distortion: ndarray

Distortion tail [alpha, lambda1, lambda2, tau_x, tau_y].

params property 

params: ndarray

Flat parameter vector in param_names order.

from_dict classmethod 

from_dict(d: dict) -> 'DSPlusModel'

Reconstruct from :meth:to_dict output.

Source code in ds_msp/models/dsplus.py 
229
230
231
232
@classmethod
def from_dict(cls, d: dict) -> "DSPlusModel":
    """Reconstruct from :meth:`to_dict` output."""
    return cls(**{k: d[k] for k in cls.param_names})

from_params classmethod 

from_params(p: ndarray) -> 'DSPlusModel'

Build from a flat vector in param_names order.

Source code in ds_msp/models/dsplus.py 
190
191
192
193
@classmethod
def from_params(cls, p: np.ndarray) -> "DSPlusModel":
    """Build from a flat vector in `param_names` order."""
    return cls(*np.asarray(p, dtype=np.float64).ravel())

initialize_from_correspondences 

initialize_from_correspondences(K_seed, rays, pixels) -> None

Seed fx,fy,cx,cy from K_seed; solve alpha linearly (same UCM closed form as UCMModel); zero the radial/tilt terms.

Source code in ds_msp/models/dsplus.py 
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
    """Seed ``fx,fy,cx,cy`` from `K_seed`; solve ``alpha`` linearly (same UCM
    closed form as `UCMModel`); zero the radial/tilt terms."""
    self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
    self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
    rays = np.asarray(rays, dtype=np.float64)
    x, y, z = rays[:, 0], rays[:, 1], rays[:, 2]
    mx = (pixels[:, 0] - self.cx) / self.fx
    my = (pixels[:, 1] - self.cy) / self.fy
    # Same UCM linear solve as ucm.py: alpha = (x - mx*z)/(mx*(1-z)).
    A = np.concatenate([mx * (1.0 - z), my * (1.0 - z)])
    b = np.concatenate([x - mx * z, y - my * z])
    denom = float(A @ A)
    self.alpha = float(np.clip((A @ b) / denom, 1e-6, 1.0 - 1e-6)) if denom > 1e-12 else 0.5
    self.lambda1 = 0.0
    self.lambda2 = 0.0
    self.tau_x = 0.0
    self.tau_y = 0.0

param_bounds classmethod 

param_bounds() -> Tuple[np.ndarray, np.ndarray]

Optimizer bounds: alpha in (0, 1), tilt/radial terms modestly bounded.

Source code in ds_msp/models/dsplus.py 
195
196
197
198
199
200
201
202
@classmethod
def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
    """Optimizer bounds: ``alpha in (0, 1)``, tilt/radial terms modestly bounded."""
    lb = np.array([1.0, 1.0, -1e5, -1e5, 1e-6, -2.0, -2.0, -0.5, -0.5],
                  dtype=np.float64)
    ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0 - 1e-6, 2.0, 2.0, 0.5, 0.5],
                  dtype=np.float64)
    return lb, ub

project 

project(P: ndarray) -> Tuple[np.ndarray, np.ndarray]

Project camera-frame points through 4 closed-form-invertible stages.

DS+ replaces Double Sphere's xi shift (the proven near-null direction for the target lens, see the module docstring) with a UCM sphere core followed by two extra closed-form stages: a Fitzgibbon division-model radial layer and a 2-axis Scheimpflug tilt homography, so the composed map stays fully closed-form invertible end to end::

pixel = K . H_tau . D_lambda . S_alpha(bearing)

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates (u, v), origin top-left.
valid (ndarray, shape(...))

True iff all 4 stages are individually valid: the UCM half-space/denominator condition, a positive division-model radial factor, and a nonzero tilt-homography denominator. See ds_msp.models.dsplus_math.dsplus_project.
References

UCM core: Geyer, C., Daniilidis, K. "A Unifying Theory for Central Panoramic Systems." ECCV 2000; Mei, C., Rives, P. ICRA 2007. Division radial layer: Fitzgibbon, A. "Simultaneous linear estimation of multiple view geometry and lens distortion." CVPR 2001. Tilt homography: cf. OpenCV CALIB_TILTED_MODEL. Staged composition and Jacobian are this repo's own extension — see ds_msp/models/dsplus_math.py and ADR-0005.

Examples:

>>> import numpy as np
>>> from ds_msp.models import DSPlusModel
>>> m = DSPlusModel.sample()
>>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
>>> np.round(uv, 2)
array([[949.18, 518.81]])
Source code in ds_msp/models/dsplus.py 
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Project camera-frame points through 4 closed-form-invertible stages.

    DS+ replaces Double Sphere's ``xi`` shift (the proven near-null
    direction for the target lens, see the module docstring) with a UCM
    sphere core followed by two extra closed-form stages: a Fitzgibbon
    division-model radial layer and a 2-axis Scheimpflug tilt homography,
    so the composed map stays fully closed-form invertible end to end::

        pixel = K . H_tau . D_lambda . S_alpha(bearing)

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates ``(u, v)``, origin top-left.
    valid : ndarray, shape (...,)
        ``True`` iff all 4 stages are individually valid: the UCM
        half-space/denominator condition, a positive division-model
        radial factor, and a nonzero tilt-homography denominator. See
        ``ds_msp.models.dsplus_math.dsplus_project``.

    References
    ----------
    UCM core: Geyer, C., Daniilidis, K. "A Unifying Theory for Central
    Panoramic Systems." ECCV 2000; Mei, C., Rives, P. ICRA 2007.
    Division radial layer: Fitzgibbon, A. "Simultaneous linear estimation
    of multiple view geometry and lens distortion." CVPR 2001.
    Tilt homography: cf. OpenCV ``CALIB_TILTED_MODEL``. Staged
    composition and Jacobian are this repo's own extension — see
    ``ds_msp/models/dsplus_math.py`` and
    [ADR-0005](https://github.com/Munna-Manoj/DS-MSP/blob/main/docs/process/architecture/decisions/ADR-0005-dsplus-eucmplus.md).

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import DSPlusModel
    >>> m = DSPlusModel.sample()
    >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
    >>> np.round(uv, 2)
    array([[949.18, 518.81]])
    """
    u, v, valid = dsplus_project(
        np.asarray(P, dtype=np.float64),
        self.fx, self.fy, self.cx, self.cy, self.alpha,
        self.lambda1, self.lambda2, self.tau_x, self.tau_y)
    return np.stack([u, v], axis=-1), valid

project_jacobian 

project_jacobian(P)

Project with analytic derivatives via the chain rule through all 4 stages.

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates, identical to project.
J_point (ndarray, shape(..., 2, 3))

d(u, v) / d(x, y, z), computed as K @ J_H @ J_D @ J_S.
J_param (ndarray, shape(..., 2, 9))

d(u, v) / d(fx, fy, cx, cy, alpha, lambda1, lambda2, tau_x, tau_y), columns in param_names order.
valid (ndarray, shape(...))

Projectability mask, identical condition to project.

Examples:

>>> import numpy as np
>>> from ds_msp.models import DSPlusModel
>>> m = DSPlusModel.sample()
>>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
>>> J_point.shape, J_param.shape
((1, 2, 3), (1, 2, 9))
Source code in ds_msp/models/dsplus.py 
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
def project_jacobian(self, P):
    """Project with analytic derivatives via the chain rule through all 4 stages.

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates, identical to `project`.
    J_point : ndarray, shape (..., 2, 3)
        ``d(u, v) / d(x, y, z)``, computed as ``K @ J_H @ J_D @ J_S``.
    J_param : ndarray, shape (..., 2, 9)
        ``d(u, v) / d(fx, fy, cx, cy, alpha, lambda1, lambda2, tau_x, tau_y)``,
        columns in `param_names` order.
    valid : ndarray, shape (...,)
        Projectability mask, identical condition to `project`.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import DSPlusModel
    >>> m = DSPlusModel.sample()
    >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
    >>> J_point.shape, J_param.shape
    ((1, 2, 3), (1, 2, 9))
    """
    u, v, J_point, J_param, valid = dsplus_project_jacobian(
        np.asarray(P, dtype=np.float64),
        self.fx, self.fy, self.cx, self.cy, self.alpha,
        self.lambda1, self.lambda2, self.tau_x, self.tau_y)
    return np.stack([u, v], axis=-1), J_point, J_param, valid

sample classmethod 

sample() -> 'DSPlusModel'

Realistic instance for contract testing (the bundled calibration).

Source code in ds_msp/models/dsplus.py 
42
43
44
45
@classmethod
def sample(cls) -> "DSPlusModel":
    """Realistic instance for contract testing (the bundled calibration)."""
    return cls(711.57, 711.24, 949.18, 518.81, 0.62, -0.10, 0.02, 0.001, -0.001)

to_dict 

to_dict() -> dict

Serialize to {"model": "dsplus", "fx": ..., ..., "tau_y": ...}.

Source code in ds_msp/models/dsplus.py 
223
224
225
226
227
def to_dict(self) -> dict:
    """Serialize to ``{"model": "dsplus", "fx": ..., ..., "tau_y": ...}``."""
    d = {"model": self.name}
    d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
    return d

unproject 

unproject(uv: ndarray) -> Tuple[np.ndarray, np.ndarray]

Unproject pixels to unit bearing rays (closed form, no iteration).

Inverts the 4-stage chain in reverse (K^-1, tilt inverse, then a quadratic/quartic radical to invert the division-model radial layer, then the UCM closed form). See ds_msp.models.dsplus_math.dsplus_unproject and ds_msp.models.dsplus_math._invert_division.

Parameters:

Name Type Description Default
uv (ndarray, shape(..., 2))

Pixel coordinates (u, v), origin top-left.

 required

Returns:

Name Type Description
rays (ndarray, shape(..., 3))

Unit-norm camera-frame bearing vectors, +Z forward.
valid (ndarray, shape(...))

True iff the tilt inverse, the division-radial root, and the UCM sphere inverse are all well-defined. Invalid rays are zeroed, never NaN.

Examples:

>>> import numpy as np
>>> from ds_msp.models import DSPlusModel
>>> m = DSPlusModel.sample()
>>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
>>> np.round(rays, 4)
array([[0., 0., 1.]])
Source code in ds_msp/models/dsplus.py 
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Unproject pixels to unit bearing rays (closed form, no iteration).

    Inverts the 4-stage chain in reverse (``K^-1``, tilt inverse, then a
    quadratic/quartic radical to invert the division-model radial layer,
    then the UCM closed form). See
    ``ds_msp.models.dsplus_math.dsplus_unproject`` and
    ``ds_msp.models.dsplus_math._invert_division``.

    Parameters
    ----------
    uv : ndarray, shape (..., 2)
        Pixel coordinates ``(u, v)``, origin top-left.

    Returns
    -------
    rays : ndarray, shape (..., 3)
        Unit-norm camera-frame bearing vectors, +Z forward.
    valid : ndarray, shape (...,)
        ``True`` iff the tilt inverse, the division-radial root, and the
        UCM sphere inverse are all well-defined. Invalid rays are
        zeroed, never NaN.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import DSPlusModel
    >>> m = DSPlusModel.sample()
    >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
    >>> np.round(rays, 4)
    array([[0., 0., 1.]])
    """
    return dsplus_unproject(
        np.asarray(uv, dtype=np.float64),
        self.fx, self.fy, self.cx, self.cy, self.alpha,
        self.lambda1, self.lambda2, self.tau_x, self.tau_y)

Extended Unified Camera Model (EUCM)

ds_msp.models.EUCMModel

Enhanced UCM (Khomutenko et al. 2016). Satisfies CameraModel.

Source code in ds_msp/models/eucm.py 
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
class EUCMModel:
    """Enhanced UCM (Khomutenko et al. 2016). Satisfies ``CameraModel``."""

    name: ClassVar[str] = "eucm"
    param_names: ClassVar[Tuple[str, ...]] = ("fx", "fy", "cx", "cy", "alpha", "beta")

    def __init__(self, fx: float, fy: float, cx: float, cy: float,
                 alpha: float, beta: float) -> None:
        self.fx = float(fx)
        self.fy = float(fy)
        self.cx = float(cx)
        self.cy = float(cy)
        self.alpha = float(alpha)
        self.beta = float(beta)

    @classmethod
    def sample(cls) -> "EUCMModel":
        """Realistic instance for contract testing (the bundled calibration)."""
        return cls(711.57, 711.24, 949.18, 518.81, 0.6, 1.1)

    @property
    def params(self) -> np.ndarray:
        """Flat parameter vector ``[fx, fy, cx, cy, alpha, beta]``."""
        return np.array([self.fx, self.fy, self.cx, self.cy, self.alpha, self.beta],
                        dtype=np.float64)

    @property
    def K(self) -> np.ndarray:
        """3x3 pinhole intrinsic matrix built from ``fx, fy, cx, cy``."""
        return np.array([[self.fx, 0.0, self.cx], [0.0, self.fy, self.cy],
                         [0.0, 0.0, 1.0]], dtype=np.float64)

    @property
    def distortion(self) -> np.ndarray:
        """Distortion tail ``[alpha, beta]`` (perspective blend, ellipse-radius weight)."""
        return np.array([self.alpha, self.beta], dtype=np.float64)

    def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Project camera-frame points via an ellipse-radius-weighted UCM sphere.

        EUCM generalizes the Unified Camera Model (`UCMModel`) with a second
        parameter ``beta`` that reweights the radial term inside the sphere
        distance, ``d = sqrt(beta*(x^2+y^2) + z^2)``, before the same
        ``alpha``-blended perspective division. This extra degree of freedom
        (decoupling the radial and axial curvature of the projection surface)
        is what lets EUCM fit distortion profiles a single-sphere UCM cannot.

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates ``(u, v)``, origin top-left.
        valid : ndarray, shape (...,)
            ``True`` iff the perspective-division denominator
            ``alpha*d + (1-alpha)*z`` is bounded away from zero. See
            ``ds_msp.models.eucm_math.eucm_project``.

        References
        ----------
        Khomutenko, B., Garcia, G., Martinet, P. "An Enhanced Unified Camera
        Model for Omnidirectional Cameras." IEEE RA-L 2016.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import EUCMModel
        >>> m = EUCMModel.sample()
        >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
        >>> np.round(uv, 2)
        array([[949.18, 518.81]])
        """
        u, v, valid = eucm_project(np.asarray(P, dtype=np.float64),
                                   self.fx, self.fy, self.cx, self.cy, self.alpha, self.beta)
        return np.stack([u, v], axis=-1), valid

    def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Unproject pixels to unit bearing rays (closed form).

        Parameters
        ----------
        uv : ndarray, shape (..., 2)
            Pixel coordinates ``(u, v)``, origin top-left.

        Returns
        -------
        rays : ndarray, shape (..., 3)
            Unit-norm camera-frame bearing vectors, +Z forward.
        valid : ndarray, shape (...,)
            ``True`` iff the pixel lies within the reachable disc and the
            perspective-division denominator is bounded away from zero.
            Invalid rays are zeroed, never NaN. See
            ``ds_msp.models.eucm_math.eucm_unproject``.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import EUCMModel
        >>> m = EUCMModel.sample()
        >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
        >>> np.round(rays, 4)
        array([[0., 0., 1.]])
        """
        return eucm_unproject(np.asarray(uv, dtype=np.float64),
                              self.fx, self.fy, self.cx, self.cy, self.alpha, self.beta)

    def project_jacobian(self, P):
        """Project with analytic derivatives (no autodiff, no finite differences).

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates, identical to `project`.
        J_point : ndarray, shape (..., 2, 3)
            ``d(u, v) / d(x, y, z)``.
        J_param : ndarray, shape (..., 2, 6)
            ``d(u, v) / d(fx, fy, cx, cy, alpha, beta)``, columns in
            `param_names` order.
        valid : ndarray, shape (...,)
            Projectability mask, identical condition to `project`.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import EUCMModel
        >>> m = EUCMModel.sample()
        >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
        >>> J_point.shape, J_param.shape
        ((1, 2, 3), (1, 2, 6))
        """
        u, v, J_point, J_param, valid = eucm_project_jacobian(
            np.asarray(P, dtype=np.float64),
            self.fx, self.fy, self.cx, self.cy, self.alpha, self.beta)
        return np.stack([u, v], axis=-1), J_point, J_param, valid

    @classmethod
    def from_params(cls, p: np.ndarray) -> "EUCMModel":
        """Build from a flat ``[fx, fy, cx, cy, alpha, beta]`` vector."""
        return cls(*np.asarray(p, dtype=np.float64).ravel())

    @classmethod
    def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
        """Optimizer bounds: ``alpha in (0, 1)``, ``beta in (0, 10]``."""
        lb = np.array([1.0, 1.0, -1e5, -1e5, 1e-6, 1e-3], dtype=np.float64)
        ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0 - 1e-6, 10.0], dtype=np.float64)
        return lb, ub

    def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
        """Seed ``fx,fy,cx,cy`` from `K_seed`; fix ``beta=1`` and solve ``alpha`` linearly
        (the UCM closed form; ``beta=1`` reduces EUCM to UCM)."""
        self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
        self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
        self.beta = 1.0
        rays = np.asarray(rays, dtype=np.float64)
        x, y, z = rays[:, 0], rays[:, 1], rays[:, 2]
        mx = (pixels[:, 0] - self.cx) / self.fx
        my = (pixels[:, 1] - self.cy) / self.fy
        # beta=1, unit rays (d=1): alpha = (x - mx*z)/(mx*(1 - z)), linear LS.
        A = np.concatenate([mx * (1.0 - z), my * (1.0 - z)])
        b = np.concatenate([x - mx * z, y - my * z])
        denom = float(A @ A)
        self.alpha = float(np.clip((A @ b) / denom, 1e-6, 1.0 - 1e-6)) if denom > 1e-12 else 0.5

    def to_dict(self) -> dict:
        """Serialize to ``{"model": "eucm", "fx": ..., ..., "beta": ...}``."""
        d = {"model": self.name}
        d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
        return d

    @classmethod
    def from_dict(cls, d: dict) -> "EUCMModel":
        """Reconstruct from :meth:`to_dict` output."""
        return cls(**{k: d[k] for k in cls.param_names})

    def __repr__(self) -> str:
        return ("EUCMModel(fx={:.3f}, fy={:.3f}, cx={:.3f}, cy={:.3f}, "
                "alpha={:.4f}, beta={:.4f})").format(
                    self.fx, self.fy, self.cx, self.cy, self.alpha, self.beta)

K property 

K: ndarray

3x3 pinhole intrinsic matrix built from fx, fy, cx, cy.

distortion property 

distortion: ndarray

Distortion tail [alpha, beta] (perspective blend, ellipse-radius weight).

params property 

params: ndarray

Flat parameter vector [fx, fy, cx, cy, alpha, beta].

from_dict classmethod 

from_dict(d: dict) -> 'EUCMModel'

Reconstruct from :meth:to_dict output.

Source code in ds_msp/models/eucm.py 
189
190
191
192
@classmethod
def from_dict(cls, d: dict) -> "EUCMModel":
    """Reconstruct from :meth:`to_dict` output."""
    return cls(**{k: d[k] for k in cls.param_names})

from_params classmethod 

from_params(p: ndarray) -> 'EUCMModel'

Build from a flat [fx, fy, cx, cy, alpha, beta] vector.

Source code in ds_msp/models/eucm.py 
155
156
157
158
@classmethod
def from_params(cls, p: np.ndarray) -> "EUCMModel":
    """Build from a flat ``[fx, fy, cx, cy, alpha, beta]`` vector."""
    return cls(*np.asarray(p, dtype=np.float64).ravel())

initialize_from_correspondences 

initialize_from_correspondences(K_seed, rays, pixels) -> None

Seed fx,fy,cx,cy from K_seed; fix beta=1 and solve alpha linearly (the UCM closed form; beta=1 reduces EUCM to UCM).

Source code in ds_msp/models/eucm.py 
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
    """Seed ``fx,fy,cx,cy`` from `K_seed`; fix ``beta=1`` and solve ``alpha`` linearly
    (the UCM closed form; ``beta=1`` reduces EUCM to UCM)."""
    self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
    self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
    self.beta = 1.0
    rays = np.asarray(rays, dtype=np.float64)
    x, y, z = rays[:, 0], rays[:, 1], rays[:, 2]
    mx = (pixels[:, 0] - self.cx) / self.fx
    my = (pixels[:, 1] - self.cy) / self.fy
    # beta=1, unit rays (d=1): alpha = (x - mx*z)/(mx*(1 - z)), linear LS.
    A = np.concatenate([mx * (1.0 - z), my * (1.0 - z)])
    b = np.concatenate([x - mx * z, y - my * z])
    denom = float(A @ A)
    self.alpha = float(np.clip((A @ b) / denom, 1e-6, 1.0 - 1e-6)) if denom > 1e-12 else 0.5

param_bounds classmethod 

param_bounds() -> Tuple[np.ndarray, np.ndarray]

Optimizer bounds: alpha in (0, 1), beta in (0, 10].

Source code in ds_msp/models/eucm.py 
160
161
162
163
164
165
@classmethod
def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
    """Optimizer bounds: ``alpha in (0, 1)``, ``beta in (0, 10]``."""
    lb = np.array([1.0, 1.0, -1e5, -1e5, 1e-6, 1e-3], dtype=np.float64)
    ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0 - 1e-6, 10.0], dtype=np.float64)
    return lb, ub

project 

project(P: ndarray) -> Tuple[np.ndarray, np.ndarray]

Project camera-frame points via an ellipse-radius-weighted UCM sphere.

EUCM generalizes the Unified Camera Model (UCMModel) with a second parameter beta that reweights the radial term inside the sphere distance, d = sqrt(beta*(x^2+y^2) + z^2), before the same alpha-blended perspective division. This extra degree of freedom (decoupling the radial and axial curvature of the projection surface) is what lets EUCM fit distortion profiles a single-sphere UCM cannot.

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates (u, v), origin top-left.
valid (ndarray, shape(...))

True iff the perspective-division denominator alpha*d + (1-alpha)*z is bounded away from zero. See ds_msp.models.eucm_math.eucm_project.
References

Khomutenko, B., Garcia, G., Martinet, P. "An Enhanced Unified Camera Model for Omnidirectional Cameras." IEEE RA-L 2016.

Examples:

>>> import numpy as np
>>> from ds_msp.models import EUCMModel
>>> m = EUCMModel.sample()
>>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
>>> np.round(uv, 2)
array([[949.18, 518.81]])
Source code in ds_msp/models/eucm.py 
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Project camera-frame points via an ellipse-radius-weighted UCM sphere.

    EUCM generalizes the Unified Camera Model (`UCMModel`) with a second
    parameter ``beta`` that reweights the radial term inside the sphere
    distance, ``d = sqrt(beta*(x^2+y^2) + z^2)``, before the same
    ``alpha``-blended perspective division. This extra degree of freedom
    (decoupling the radial and axial curvature of the projection surface)
    is what lets EUCM fit distortion profiles a single-sphere UCM cannot.

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates ``(u, v)``, origin top-left.
    valid : ndarray, shape (...,)
        ``True`` iff the perspective-division denominator
        ``alpha*d + (1-alpha)*z`` is bounded away from zero. See
        ``ds_msp.models.eucm_math.eucm_project``.

    References
    ----------
    Khomutenko, B., Garcia, G., Martinet, P. "An Enhanced Unified Camera
    Model for Omnidirectional Cameras." IEEE RA-L 2016.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import EUCMModel
    >>> m = EUCMModel.sample()
    >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
    >>> np.round(uv, 2)
    array([[949.18, 518.81]])
    """
    u, v, valid = eucm_project(np.asarray(P, dtype=np.float64),
                               self.fx, self.fy, self.cx, self.cy, self.alpha, self.beta)
    return np.stack([u, v], axis=-1), valid

project_jacobian 

project_jacobian(P)

Project with analytic derivatives (no autodiff, no finite differences).

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates, identical to project.
J_point (ndarray, shape(..., 2, 3))

d(u, v) / d(x, y, z).
J_param (ndarray, shape(..., 2, 6))

d(u, v) / d(fx, fy, cx, cy, alpha, beta), columns in param_names order.
valid (ndarray, shape(...))

Projectability mask, identical condition to project.

Examples:

>>> import numpy as np
>>> from ds_msp.models import EUCMModel
>>> m = EUCMModel.sample()
>>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
>>> J_point.shape, J_param.shape
((1, 2, 3), (1, 2, 6))
Source code in ds_msp/models/eucm.py 
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
def project_jacobian(self, P):
    """Project with analytic derivatives (no autodiff, no finite differences).

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates, identical to `project`.
    J_point : ndarray, shape (..., 2, 3)
        ``d(u, v) / d(x, y, z)``.
    J_param : ndarray, shape (..., 2, 6)
        ``d(u, v) / d(fx, fy, cx, cy, alpha, beta)``, columns in
        `param_names` order.
    valid : ndarray, shape (...,)
        Projectability mask, identical condition to `project`.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import EUCMModel
    >>> m = EUCMModel.sample()
    >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
    >>> J_point.shape, J_param.shape
    ((1, 2, 3), (1, 2, 6))
    """
    u, v, J_point, J_param, valid = eucm_project_jacobian(
        np.asarray(P, dtype=np.float64),
        self.fx, self.fy, self.cx, self.cy, self.alpha, self.beta)
    return np.stack([u, v], axis=-1), J_point, J_param, valid

sample classmethod 

sample() -> 'EUCMModel'

Realistic instance for contract testing (the bundled calibration).

Source code in ds_msp/models/eucm.py 
27
28
29
30
@classmethod
def sample(cls) -> "EUCMModel":
    """Realistic instance for contract testing (the bundled calibration)."""
    return cls(711.57, 711.24, 949.18, 518.81, 0.6, 1.1)

to_dict 

to_dict() -> dict

Serialize to {"model": "eucm", "fx": ..., ..., "beta": ...}.

Source code in ds_msp/models/eucm.py 
183
184
185
186
187
def to_dict(self) -> dict:
    """Serialize to ``{"model": "eucm", "fx": ..., ..., "beta": ...}``."""
    d = {"model": self.name}
    d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
    return d

unproject 

unproject(uv: ndarray) -> Tuple[np.ndarray, np.ndarray]

Unproject pixels to unit bearing rays (closed form).

Parameters:

Name Type Description Default
uv (ndarray, shape(..., 2))

Pixel coordinates (u, v), origin top-left.

 required

Returns:

Name Type Description
rays (ndarray, shape(..., 3))

Unit-norm camera-frame bearing vectors, +Z forward.
valid (ndarray, shape(...))

True iff the pixel lies within the reachable disc and the perspective-division denominator is bounded away from zero. Invalid rays are zeroed, never NaN. See ds_msp.models.eucm_math.eucm_unproject.

Examples:

>>> import numpy as np
>>> from ds_msp.models import EUCMModel
>>> m = EUCMModel.sample()
>>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
>>> np.round(rays, 4)
array([[0., 0., 1.]])
Source code in ds_msp/models/eucm.py 
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Unproject pixels to unit bearing rays (closed form).

    Parameters
    ----------
    uv : ndarray, shape (..., 2)
        Pixel coordinates ``(u, v)``, origin top-left.

    Returns
    -------
    rays : ndarray, shape (..., 3)
        Unit-norm camera-frame bearing vectors, +Z forward.
    valid : ndarray, shape (...,)
        ``True`` iff the pixel lies within the reachable disc and the
        perspective-division denominator is bounded away from zero.
        Invalid rays are zeroed, never NaN. See
        ``ds_msp.models.eucm_math.eucm_unproject``.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import EUCMModel
    >>> m = EUCMModel.sample()
    >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
    >>> np.round(rays, 4)
    array([[0., 0., 1.]])
    """
    return eucm_unproject(np.asarray(uv, dtype=np.float64),
                          self.fx, self.fy, self.cx, self.cy, self.alpha, self.beta)

Kannala-Brandt (KB)

ds_msp.models.KannalaBrandtModel

Kannala-Brandt / equidistant fisheye. Satisfies CameraModel.

K and distortion ([k1,k2,k3,k4]) plug directly into cv2.fisheye.

Source code in ds_msp/models/kb.py 
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
class KannalaBrandtModel:
    """Kannala-Brandt / equidistant fisheye. Satisfies ``CameraModel``.

    ``K`` and ``distortion`` ([k1,k2,k3,k4]) plug directly into ``cv2.fisheye``.
    """

    name: ClassVar[str] = "kb"
    param_names: ClassVar[Tuple[str, ...]] = (
        "fx", "fy", "cx", "cy", "k1", "k2", "k3", "k4")

    def __init__(self, fx, fy, cx, cy, k1=0.0, k2=0.0, k3=0.0, k4=0.0) -> None:
        self.fx = float(fx)
        self.fy = float(fy)
        self.cx = float(cx)
        self.cy = float(cy)
        self.k1 = float(k1)
        self.k2 = float(k2)
        self.k3 = float(k3)
        self.k4 = float(k4)

    @classmethod
    def sample(cls) -> "KannalaBrandtModel":
        """Realistic instance for contract testing (a narrower-FOV KB lens)."""
        return cls(320.0, 321.0, 320.0, 240.0, 0.05, 0.01, -0.002, 0.0008)

    @property
    def params(self) -> np.ndarray:
        """Flat parameter vector ``[fx, fy, cx, cy, k1, k2, k3, k4]``."""
        return np.array([self.fx, self.fy, self.cx, self.cy,
                         self.k1, self.k2, self.k3, self.k4], dtype=np.float64)

    @property
    def K(self) -> np.ndarray:
        """3x3 pinhole intrinsic matrix built from ``fx, fy, cx, cy``."""
        return np.array([[self.fx, 0.0, self.cx], [0.0, self.fy, self.cy],
                         [0.0, 0.0, 1.0]], dtype=np.float64)

    @property
    def distortion(self) -> np.ndarray:
        """Odd-power angle-polynomial coefficients ``[k1, k2, k3, k4]`` (OpenCV order)."""
        return np.array([self.k1, self.k2, self.k3, self.k4], dtype=np.float64)

    def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Project camera-frame points via the equidistant angle-polynomial map.

        Kannala-Brandt is distinct in operating on the incidence angle rather
        than a sphere: it computes ``theta = atan2(r, z)`` (the angle off the
        optical axis) and maps it through an odd-power polynomial
        ``theta_d = theta + k1*theta^3 + k2*theta^5 + k3*theta^7 +
        k4*theta^9``, then scales the direction ``(x, y)/r`` by
        ``theta_d / r``. With ``k1..k4 = 0`` this reduces to the ideal
        equidistant fisheye ``r = f*theta``.

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates ``(u, v)``, origin top-left.
        valid : ndarray, shape (...,)
            ``True`` iff both pixel coordinates are finite. See
            ``ds_msp.models.kb_math.kb_project``.

        References
        ----------
        Kannala, J., Brandt, S. S. "A Generic Camera Model and Calibration
        Method for Conventional, Wide-Angle, and Fish-Eye Lenses." IEEE
        TPAMI 2006. ``K`` and `distortion` plug directly into
        ``cv2.fisheye.projectPoints``.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import KannalaBrandtModel
        >>> m = KannalaBrandtModel.sample()
        >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
        >>> np.round(uv, 2)
        array([[320., 240.]])
        """
        u, v, valid = kb_project(np.asarray(P, dtype=np.float64),
                                 self.fx, self.fy, self.cx, self.cy,
                                 self.k1, self.k2, self.k3, self.k4)
        return np.stack([u, v], axis=-1), valid

    def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Unproject pixels to unit bearing rays via Newton-Raphson.

        Inverts ``theta_d(theta) = ru`` for ``theta`` with 10 fixed Newton
        iterations (no closed form exists for a degree-9 polynomial), then
        builds the ray from ``(sin(theta), cos(theta))`` and the pixel's
        angular direction. See ``ds_msp.models.kb_math.kb_unproject``.

        Parameters
        ----------
        uv : ndarray, shape (..., 2)
            Pixel coordinates ``(u, v)``, origin top-left.

        Returns
        -------
        rays : ndarray, shape (..., 3)
            Unit-norm camera-frame bearing vectors, +Z forward.
        valid : ndarray, shape (...,)
            ``True`` iff the Newton residual ``|theta_d(theta) - ru| < 1e-6``
            and ``theta <= pi``. Invalid rays are zeroed, never NaN.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import KannalaBrandtModel
        >>> m = KannalaBrandtModel.sample()
        >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
        >>> np.round(rays, 4)
        array([[0., 0., 1.]])
        """
        return kb_unproject(np.asarray(uv, dtype=np.float64),
                            self.fx, self.fy, self.cx, self.cy,
                            self.k1, self.k2, self.k3, self.k4)

    def project_jacobian(self, P):
        """Project with analytic derivatives (no autodiff, no finite differences).

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates, identical to `project`.
        J_point : ndarray, shape (..., 2, 3)
            ``d(u, v) / d(x, y, z)``.
        J_param : ndarray, shape (..., 2, 8)
            ``d(u, v) / d(fx, fy, cx, cy, k1, k2, k3, k4)``, columns in
            `param_names` order (``d theta_d / d k_i = theta^(2i+1)``).
        valid : ndarray, shape (...,)
            Projectability mask, identical condition to `project`.

        References
        ----------
        Kannala, J., Brandt, S. S. IEEE TPAMI 2006 (closed-form Jacobian
        derived from the forward map; verified here by finite-difference
        check, relative error <= 1e-6, see ``pytest -m jac``).

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import KannalaBrandtModel
        >>> m = KannalaBrandtModel.sample()
        >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
        >>> J_point.shape, J_param.shape
        ((1, 2, 3), (1, 2, 8))
        """
        u, v, J_point, J_param, valid = kb_project_jacobian(
            np.asarray(P, dtype=np.float64),
            self.fx, self.fy, self.cx, self.cy,
            self.k1, self.k2, self.k3, self.k4)
        return np.stack([u, v], axis=-1), J_point, J_param, valid

    @classmethod
    def from_params(cls, p: np.ndarray) -> "KannalaBrandtModel":
        """Build from a flat ``[fx, fy, cx, cy, k1, k2, k3, k4]`` vector."""
        return cls(*np.asarray(p, dtype=np.float64).ravel())

    @classmethod
    def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
        """Optimizer bounds: each ``k_i in [-1, 1]``, focal/center wide-open."""
        lb = np.array([1.0, 1.0, -1e5, -1e5, -1.0, -1.0, -1.0, -1.0], dtype=np.float64)
        ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0, 1.0, 1.0, 1.0], dtype=np.float64)
        return lb, ub

    def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
        """Seed ``fx,fy,cx,cy`` from `K_seed`; solve ``k1..k4`` by linear least squares
        (``ru - theta`` is linear in the odd powers of ``theta``)."""
        self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
        self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
        rays = np.asarray(rays, dtype=np.float64)
        theta = np.arctan2(np.sqrt(rays[:, 0]**2 + rays[:, 1]**2), rays[:, 2])
        mx = (pixels[:, 0] - self.cx) / self.fx
        my = (pixels[:, 1] - self.cy) / self.fy
        ru = np.sqrt(mx*mx + my*my)
        # ru = theta + k1 th^3 + k2 th^5 + k3 th^7 + k4 th^9 -> linear in k.
        A = np.stack([theta**3, theta**5, theta**7, theta**9], axis=1)
        b = ru - theta
        coeffs, *_ = np.linalg.lstsq(A, b, rcond=None)
        self.k1, self.k2, self.k3, self.k4 = (float(c) for c in coeffs)

    def to_dict(self) -> dict:
        """Serialize to ``{"model": "kb", "fx": ..., ..., "k4": ...}``."""
        d = {"model": self.name}
        d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
        return d

    @classmethod
    def from_dict(cls, d: dict) -> "KannalaBrandtModel":
        """Reconstruct from :meth:`to_dict` output."""
        return cls(**{k: d[k] for k in cls.param_names})

    def __repr__(self) -> str:
        return ("KannalaBrandtModel(fx={:.3f}, fy={:.3f}, cx={:.3f}, cy={:.3f}, "
                "k=[{:.5f}, {:.5f}, {:.5f}, {:.5f}])").format(
                    self.fx, self.fy, self.cx, self.cy,
                    self.k1, self.k2, self.k3, self.k4)

K property 

K: ndarray

3x3 pinhole intrinsic matrix built from fx, fy, cx, cy.

distortion property 

distortion: ndarray

Odd-power angle-polynomial coefficients [k1, k2, k3, k4] (OpenCV order).

params property 

params: ndarray

Flat parameter vector [fx, fy, cx, cy, k1, k2, k3, k4].

from_dict classmethod 

from_dict(d: dict) -> 'KannalaBrandtModel'

Reconstruct from :meth:to_dict output.

Source code in ds_msp/models/kb.py 
208
209
210
211
@classmethod
def from_dict(cls, d: dict) -> "KannalaBrandtModel":
    """Reconstruct from :meth:`to_dict` output."""
    return cls(**{k: d[k] for k in cls.param_names})

from_params classmethod 

from_params(p: ndarray) -> 'KannalaBrandtModel'

Build from a flat [fx, fy, cx, cy, k1, k2, k3, k4] vector.

Source code in ds_msp/models/kb.py 
174
175
176
177
@classmethod
def from_params(cls, p: np.ndarray) -> "KannalaBrandtModel":
    """Build from a flat ``[fx, fy, cx, cy, k1, k2, k3, k4]`` vector."""
    return cls(*np.asarray(p, dtype=np.float64).ravel())

initialize_from_correspondences 

initialize_from_correspondences(K_seed, rays, pixels) -> None

Seed fx,fy,cx,cy from K_seed; solve k1..k4 by linear least squares (ru - theta is linear in the odd powers of theta).

Source code in ds_msp/models/kb.py 
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
    """Seed ``fx,fy,cx,cy`` from `K_seed`; solve ``k1..k4`` by linear least squares
    (``ru - theta`` is linear in the odd powers of ``theta``)."""
    self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
    self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
    rays = np.asarray(rays, dtype=np.float64)
    theta = np.arctan2(np.sqrt(rays[:, 0]**2 + rays[:, 1]**2), rays[:, 2])
    mx = (pixels[:, 0] - self.cx) / self.fx
    my = (pixels[:, 1] - self.cy) / self.fy
    ru = np.sqrt(mx*mx + my*my)
    # ru = theta + k1 th^3 + k2 th^5 + k3 th^7 + k4 th^9 -> linear in k.
    A = np.stack([theta**3, theta**5, theta**7, theta**9], axis=1)
    b = ru - theta
    coeffs, *_ = np.linalg.lstsq(A, b, rcond=None)
    self.k1, self.k2, self.k3, self.k4 = (float(c) for c in coeffs)

param_bounds classmethod 

param_bounds() -> Tuple[np.ndarray, np.ndarray]

Optimizer bounds: each k_i in [-1, 1], focal/center wide-open.

Source code in ds_msp/models/kb.py 
179
180
181
182
183
184
@classmethod
def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
    """Optimizer bounds: each ``k_i in [-1, 1]``, focal/center wide-open."""
    lb = np.array([1.0, 1.0, -1e5, -1e5, -1.0, -1.0, -1.0, -1.0], dtype=np.float64)
    ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0, 1.0, 1.0, 1.0], dtype=np.float64)
    return lb, ub

project 

project(P: ndarray) -> Tuple[np.ndarray, np.ndarray]

Project camera-frame points via the equidistant angle-polynomial map.

Kannala-Brandt is distinct in operating on the incidence angle rather than a sphere: it computes theta = atan2(r, z) (the angle off the optical axis) and maps it through an odd-power polynomial theta_d = theta + k1*theta^3 + k2*theta^5 + k3*theta^7 + k4*theta^9, then scales the direction (x, y)/r by theta_d / r. With k1..k4 = 0 this reduces to the ideal equidistant fisheye r = f*theta.

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates (u, v), origin top-left.
valid (ndarray, shape(...))

True iff both pixel coordinates are finite. See ds_msp.models.kb_math.kb_project.
References

Kannala, J., Brandt, S. S. "A Generic Camera Model and Calibration Method for Conventional, Wide-Angle, and Fish-Eye Lenses." IEEE TPAMI 2006. K and distortion plug directly into cv2.fisheye.projectPoints.

Examples:

>>> import numpy as np
>>> from ds_msp.models import KannalaBrandtModel
>>> m = KannalaBrandtModel.sample()
>>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
>>> np.round(uv, 2)
array([[320., 240.]])
Source code in ds_msp/models/kb.py 
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Project camera-frame points via the equidistant angle-polynomial map.

    Kannala-Brandt is distinct in operating on the incidence angle rather
    than a sphere: it computes ``theta = atan2(r, z)`` (the angle off the
    optical axis) and maps it through an odd-power polynomial
    ``theta_d = theta + k1*theta^3 + k2*theta^5 + k3*theta^7 +
    k4*theta^9``, then scales the direction ``(x, y)/r`` by
    ``theta_d / r``. With ``k1..k4 = 0`` this reduces to the ideal
    equidistant fisheye ``r = f*theta``.

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates ``(u, v)``, origin top-left.
    valid : ndarray, shape (...,)
        ``True`` iff both pixel coordinates are finite. See
        ``ds_msp.models.kb_math.kb_project``.

    References
    ----------
    Kannala, J., Brandt, S. S. "A Generic Camera Model and Calibration
    Method for Conventional, Wide-Angle, and Fish-Eye Lenses." IEEE
    TPAMI 2006. ``K`` and `distortion` plug directly into
    ``cv2.fisheye.projectPoints``.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import KannalaBrandtModel
    >>> m = KannalaBrandtModel.sample()
    >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
    >>> np.round(uv, 2)
    array([[320., 240.]])
    """
    u, v, valid = kb_project(np.asarray(P, dtype=np.float64),
                             self.fx, self.fy, self.cx, self.cy,
                             self.k1, self.k2, self.k3, self.k4)
    return np.stack([u, v], axis=-1), valid

project_jacobian 

project_jacobian(P)

Project with analytic derivatives (no autodiff, no finite differences).

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates, identical to project.
J_point (ndarray, shape(..., 2, 3))

d(u, v) / d(x, y, z).
J_param (ndarray, shape(..., 2, 8))

d(u, v) / d(fx, fy, cx, cy, k1, k2, k3, k4), columns in param_names order (d theta_d / d k_i = theta^(2i+1)).
valid (ndarray, shape(...))

Projectability mask, identical condition to project.
References

Kannala, J., Brandt, S. S. IEEE TPAMI 2006 (closed-form Jacobian derived from the forward map; verified here by finite-difference check, relative error <= 1e-6, see pytest -m jac).

Examples:

>>> import numpy as np
>>> from ds_msp.models import KannalaBrandtModel
>>> m = KannalaBrandtModel.sample()
>>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
>>> J_point.shape, J_param.shape
((1, 2, 3), (1, 2, 8))
Source code in ds_msp/models/kb.py 
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
def project_jacobian(self, P):
    """Project with analytic derivatives (no autodiff, no finite differences).

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates, identical to `project`.
    J_point : ndarray, shape (..., 2, 3)
        ``d(u, v) / d(x, y, z)``.
    J_param : ndarray, shape (..., 2, 8)
        ``d(u, v) / d(fx, fy, cx, cy, k1, k2, k3, k4)``, columns in
        `param_names` order (``d theta_d / d k_i = theta^(2i+1)``).
    valid : ndarray, shape (...,)
        Projectability mask, identical condition to `project`.

    References
    ----------
    Kannala, J., Brandt, S. S. IEEE TPAMI 2006 (closed-form Jacobian
    derived from the forward map; verified here by finite-difference
    check, relative error <= 1e-6, see ``pytest -m jac``).

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import KannalaBrandtModel
    >>> m = KannalaBrandtModel.sample()
    >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
    >>> J_point.shape, J_param.shape
    ((1, 2, 3), (1, 2, 8))
    """
    u, v, J_point, J_param, valid = kb_project_jacobian(
        np.asarray(P, dtype=np.float64),
        self.fx, self.fy, self.cx, self.cy,
        self.k1, self.k2, self.k3, self.k4)
    return np.stack([u, v], axis=-1), J_point, J_param, valid

sample classmethod 

sample() -> 'KannalaBrandtModel'

Realistic instance for contract testing (a narrower-FOV KB lens).

Source code in ds_msp/models/kb.py 
32
33
34
35
@classmethod
def sample(cls) -> "KannalaBrandtModel":
    """Realistic instance for contract testing (a narrower-FOV KB lens)."""
    return cls(320.0, 321.0, 320.0, 240.0, 0.05, 0.01, -0.002, 0.0008)

to_dict 

to_dict() -> dict

Serialize to {"model": "kb", "fx": ..., ..., "k4": ...}.

Source code in ds_msp/models/kb.py 
202
203
204
205
206
def to_dict(self) -> dict:
    """Serialize to ``{"model": "kb", "fx": ..., ..., "k4": ...}``."""
    d = {"model": self.name}
    d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
    return d

unproject 

unproject(uv: ndarray) -> Tuple[np.ndarray, np.ndarray]

Unproject pixels to unit bearing rays via Newton-Raphson.

Inverts theta_d(theta) = ru for theta with 10 fixed Newton iterations (no closed form exists for a degree-9 polynomial), then builds the ray from (sin(theta), cos(theta)) and the pixel's angular direction. See ds_msp.models.kb_math.kb_unproject.

Parameters:

Name Type Description Default
uv (ndarray, shape(..., 2))

Pixel coordinates (u, v), origin top-left.

 required

Returns:

Name Type Description
rays (ndarray, shape(..., 3))

Unit-norm camera-frame bearing vectors, +Z forward.
valid (ndarray, shape(...))

True iff the Newton residual |theta_d(theta) - ru| < 1e-6 and theta <= pi. Invalid rays are zeroed, never NaN.

Examples:

>>> import numpy as np
>>> from ds_msp.models import KannalaBrandtModel
>>> m = KannalaBrandtModel.sample()
>>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
>>> np.round(rays, 4)
array([[0., 0., 1.]])
Source code in ds_msp/models/kb.py 
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Unproject pixels to unit bearing rays via Newton-Raphson.

    Inverts ``theta_d(theta) = ru`` for ``theta`` with 10 fixed Newton
    iterations (no closed form exists for a degree-9 polynomial), then
    builds the ray from ``(sin(theta), cos(theta))`` and the pixel's
    angular direction. See ``ds_msp.models.kb_math.kb_unproject``.

    Parameters
    ----------
    uv : ndarray, shape (..., 2)
        Pixel coordinates ``(u, v)``, origin top-left.

    Returns
    -------
    rays : ndarray, shape (..., 3)
        Unit-norm camera-frame bearing vectors, +Z forward.
    valid : ndarray, shape (...,)
        ``True`` iff the Newton residual ``|theta_d(theta) - ru| < 1e-6``
        and ``theta <= pi``. Invalid rays are zeroed, never NaN.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import KannalaBrandtModel
    >>> m = KannalaBrandtModel.sample()
    >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
    >>> np.round(rays, 4)
    array([[0., 0., 1.]])
    """
    return kb_unproject(np.asarray(uv, dtype=np.float64),
                        self.fx, self.fy, self.cx, self.cy,
                        self.k1, self.k2, self.k3, self.k4)

OCam / Scaramuzza

ds_msp.models.OCamModel

Scaramuzza polynomial omni-camera. Satisfies CameraModel.

Parameters: centre (cx, cy), affine (c, d, e), world polynomial (a0..a4). Has no native fx/fy; K exposes a pinhole-equivalent focal ~ |a0|.

Source code in ds_msp/models/ocam.py 
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
class OCamModel:
    """Scaramuzza polynomial omni-camera. Satisfies ``CameraModel``.

    Parameters: centre (cx, cy), affine (c, d, e), world polynomial (a0..a4).
    Has no native fx/fy; ``K`` exposes a pinhole-equivalent focal ~ |a0|.
    """

    name: ClassVar[str] = "ocam"
    param_names: ClassVar[Tuple[str, ...]] = (
        "cx", "cy", "c", "d", "e", "a0", "a1", "a2", "a3", "a4")

    def __init__(self, cx, cy, c=1.0, d=0.0, e=0.0,
                 a0=-200.0, a1=0.0, a2=0.0, a3=0.0, a4=0.0) -> None:
        self.cx, self.cy = float(cx), float(cy)
        self.c, self.d, self.e = float(c), float(d), float(e)
        self.a0, self.a1, self.a2, self.a3, self.a4 = map(float, (a0, a1, a2, a3, a4))

    @classmethod
    def sample(cls) -> "OCamModel":
        """Realistic instance for contract testing (a synthetic omni-camera).

        Polynomial argument is normalized by ``R_REF=100``, so ``a2..a4``
        are O(1).
        """
        return cls(320.0, 240.0, 1.0, 0.0, 0.0, -220.0, 0.0, 6.0, -1.5, 0.25)

    def _p(self):
        return (self.cx, self.cy, self.c, self.d, self.e,
                self.a0, self.a1, self.a2, self.a3, self.a4)

    @property
    def params(self) -> np.ndarray:
        """Flat parameter vector ``[cx, cy, c, d, e, a0, a1, a2, a3, a4]``."""
        return np.array(self._p(), dtype=np.float64)

    @property
    def K(self) -> np.ndarray:
        """Pinhole-equivalent 3x3 intrinsic matrix; focal is ``|a0|`` (OCam has no native fx/fy)."""
        f = abs(self.a0)
        return np.array([[f, 0.0, self.cx], [0.0, f, self.cy], [0.0, 0.0, 1.0]],
                        dtype=np.float64)

    @property
    def distortion(self) -> np.ndarray:
        """Distortion tail ``[c, d, e, a0, a1, a2, a3, a4]`` (affine + world polynomial)."""
        return np.array([self.c, self.d, self.e,
                         self.a0, self.a1, self.a2, self.a3, self.a4], dtype=np.float64)

    def project(self, P):
        """Project camera-frame points via the Scaramuzza world polynomial.

        OCam has no closed-form projection: instead of a distortion formula
        applied to a normalized ray, it solves for the radial image distance
        ``rho`` satisfying ``w(rho) + m*rho = 0`` by Newton iteration, where
        ``w`` is a degree-4 "world" polynomial in ``rho / R_REF`` and
        ``m = Z / sqrt(X^2 + Y^2)``. The 2x2 affine stretch ``[[c, d], [e,
        1]]`` is then applied to model sensor non-squareness / skew before
        translation by the principal point. This general polynomial radial
        profile (rather than a fixed functional form) is what lets OCam fit
        distortion curves other closed-form models cannot.

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates ``(u, v)``, origin top-left.
            On-axis points (``sqrt(X^2+Y^2) < 1e-12``) map to the principal
            point.
        valid : ndarray, shape (...,)
            ``True`` iff the Newton solve for ``rho`` converged
            (``|w(rho) + m*rho| < 1e-6``) or the point is on-axis. See
            ``ds_msp.models.ocam_math.ocam_project``.

        References
        ----------
        Scaramuzza, D., Martinelli, A., Siegwart, R. "A Toolbox for Easily
        Calibrating Omnidirectional Cameras." IROS 2006.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import OCamModel
        >>> m = OCamModel.sample()
        >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
        >>> np.round(uv, 2)
        array([[320., 240.]])
        """
        u, v, valid = ocam_project(np.asarray(P, dtype=np.float64), *self._p())
        return np.stack([u, v], axis=-1), valid

    def unproject(self, uv):
        """Unproject pixels to unit bearing rays (direct polynomial evaluation).

        Inverts the affine stretch linearly, then evaluates the world
        polynomial ``w(rho)`` directly (no root-finding needed in this
        direction) to form ``ray = normalize([x, y, -w(rho)])``. See
        ``ds_msp.models.ocam_math.ocam_unproject``.

        Parameters
        ----------
        uv : ndarray, shape (..., 2)
            Pixel coordinates ``(u, v)``, origin top-left.

        Returns
        -------
        rays : ndarray, shape (..., 3)
            Unit-norm camera-frame bearing vectors, +Z forward.
        valid : ndarray, shape (...,)
            ``True`` iff the resulting ray is finite. Invalid rays are
            zeroed, never NaN.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import OCamModel
        >>> m = OCamModel.sample()
        >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
        >>> np.round(rays, 4)
        array([[0., 0., 1.]])
        """
        return ocam_unproject(np.asarray(uv, dtype=np.float64), *self._p())

    def project_jacobian(self, P):
        """Project with analytic derivatives via the implicit function theorem.

        ``rho`` is defined implicitly by ``w(rho) + m*rho = 0``; its
        derivatives w.r.t. both the 3D point and the polynomial coefficients
        are obtained without differentiating through the Newton loop, via
        the implicit function theorem (``drho = -dF / (dF/drho)``).

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates, identical to `project`.
        J_point : ndarray, shape (..., 2, 3)
            ``d(u, v) / d(X, Y, Z)``.
        J_param : ndarray, shape (..., 2, 10)
            ``d(u, v) / d(cx, cy, c, d, e, a0, a1, a2, a3, a4)``, columns in
            `param_names` order.
        valid : ndarray, shape (...,)
            Projectability mask, identical condition to `project`.

        References
        ----------
        Scaramuzza, D., Martinelli, A., Siegwart, R. IROS 2006 (Jacobian is
        this repo's own derivation via the implicit function theorem;
        verified by finite-difference check, relative error <= 1e-6, see
        ``pytest -m jac``).

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import OCamModel
        >>> m = OCamModel.sample()
        >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
        >>> J_point.shape, J_param.shape
        ((1, 2, 3), (1, 2, 10))
        """
        u, v, J_point, J_param, valid = ocam_project_jacobian(
            np.asarray(P, dtype=np.float64), *self._p())
        return np.stack([u, v], axis=-1), J_point, J_param, valid

    @classmethod
    def from_params(cls, p):
        """Build from a flat ``[cx, cy, c, d, e, a0, a1, a2, a3, a4]`` vector."""
        return cls(*np.asarray(p, dtype=np.float64).ravel())

    @classmethod
    def param_bounds(cls):
        """Optimizer bounds: ``c in (0, 10]``, ``d, e in [-1, 1]``, ``a0 < 0`` (forward-facing focal sign)."""
        lb = np.array([-1e5, -1e5, 0.1, -1.0, -1.0, -1e5, -1e3, -1e3, -1e3, -1e3],
                      dtype=np.float64)
        ub = np.array([1e5, 1e5, 10.0, 1.0, 1.0, -1.0, 1e3, 1e3, 1e3, 1e3],
                      dtype=np.float64)
        return lb, ub

    def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
        """Seed ``cx,cy`` from `K_seed`; fix the affine to identity and solve
        ``a0..a4`` by linear least squares against the observed world-polynomial values."""
        self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
        self.c, self.d, self.e = 1.0, 0.0, 0.0
        rays = np.asarray(rays, dtype=np.float64)
        x = pixels[:, 0] - self.cx
        y = pixels[:, 1] - self.cy
        rho = np.sqrt(x*x + y*y)
        rxy = np.maximum(np.sqrt(rays[:, 0]**2 + rays[:, 1]**2), 1e-9)
        # ray = [x, y, -w]/n  =>  w(rho) = -rho * rz / |rxy|
        w = -rho * rays[:, 2] / rxy
        rn = rho / R_REF
        A = np.stack([np.ones_like(rn), rn, rn**2, rn**3, rn**4], axis=1)
        coeffs, *_ = np.linalg.lstsq(A, w, rcond=None)
        self.a0, self.a1, self.a2, self.a3, self.a4 = (float(v) for v in coeffs)

    def to_dict(self) -> dict:
        """Serialize to ``{"model": "ocam", "cx": ..., ..., "a4": ...}``."""
        d = {"model": self.name}
        d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
        return d

    @classmethod
    def from_dict(cls, d: dict) -> "OCamModel":
        """Reconstruct from :meth:`to_dict` output."""
        return cls(**{k: d[k] for k in cls.param_names})

    def __repr__(self) -> str:
        return ("OCamModel(cx={:.2f}, cy={:.2f}, affine=[{:.3f},{:.3f},{:.3f}], "
                "a=[{:.4g},{:.4g},{:.4g},{:.4g},{:.4g}])").format(
                    self.cx, self.cy, self.c, self.d, self.e,
                    self.a0, self.a1, self.a2, self.a3, self.a4)

K property 

K: ndarray

Pinhole-equivalent 3x3 intrinsic matrix; focal is |a0| (OCam has no native fx/fy).

distortion property 

distortion: ndarray

Distortion tail [c, d, e, a0, a1, a2, a3, a4] (affine + world polynomial).

params property 

params: ndarray

Flat parameter vector [cx, cy, c, d, e, a0, a1, a2, a3, a4].

from_dict classmethod 

from_dict(d: dict) -> 'OCamModel'

Reconstruct from :meth:to_dict output.

Source code in ds_msp/models/ocam.py 
220
221
222
223
@classmethod
def from_dict(cls, d: dict) -> "OCamModel":
    """Reconstruct from :meth:`to_dict` output."""
    return cls(**{k: d[k] for k in cls.param_names})

from_params classmethod 

from_params(p)

Build from a flat [cx, cy, c, d, e, a0, a1, a2, a3, a4] vector.

Source code in ds_msp/models/ocam.py 
183
184
185
186
@classmethod
def from_params(cls, p):
    """Build from a flat ``[cx, cy, c, d, e, a0, a1, a2, a3, a4]`` vector."""
    return cls(*np.asarray(p, dtype=np.float64).ravel())

initialize_from_correspondences 

initialize_from_correspondences(K_seed, rays, pixels) -> None

Seed cx,cy from K_seed; fix the affine to identity and solve a0..a4 by linear least squares against the observed world-polynomial values.

Source code in ds_msp/models/ocam.py 
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
    """Seed ``cx,cy`` from `K_seed`; fix the affine to identity and solve
    ``a0..a4`` by linear least squares against the observed world-polynomial values."""
    self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
    self.c, self.d, self.e = 1.0, 0.0, 0.0
    rays = np.asarray(rays, dtype=np.float64)
    x = pixels[:, 0] - self.cx
    y = pixels[:, 1] - self.cy
    rho = np.sqrt(x*x + y*y)
    rxy = np.maximum(np.sqrt(rays[:, 0]**2 + rays[:, 1]**2), 1e-9)
    # ray = [x, y, -w]/n  =>  w(rho) = -rho * rz / |rxy|
    w = -rho * rays[:, 2] / rxy
    rn = rho / R_REF
    A = np.stack([np.ones_like(rn), rn, rn**2, rn**3, rn**4], axis=1)
    coeffs, *_ = np.linalg.lstsq(A, w, rcond=None)
    self.a0, self.a1, self.a2, self.a3, self.a4 = (float(v) for v in coeffs)

param_bounds classmethod 

param_bounds()

Optimizer bounds: c in (0, 10], d, e in [-1, 1], a0 < 0 (forward-facing focal sign).

Source code in ds_msp/models/ocam.py 
188
189
190
191
192
193
194
195
@classmethod
def param_bounds(cls):
    """Optimizer bounds: ``c in (0, 10]``, ``d, e in [-1, 1]``, ``a0 < 0`` (forward-facing focal sign)."""
    lb = np.array([-1e5, -1e5, 0.1, -1.0, -1.0, -1e5, -1e3, -1e3, -1e3, -1e3],
                  dtype=np.float64)
    ub = np.array([1e5, 1e5, 10.0, 1.0, 1.0, -1.0, 1e3, 1e3, 1e3, 1e3],
                  dtype=np.float64)
    return lb, ub

project 

project(P)

Project camera-frame points via the Scaramuzza world polynomial.

OCam has no closed-form projection: instead of a distortion formula applied to a normalized ray, it solves for the radial image distance rho satisfying w(rho) + m*rho = 0 by Newton iteration, where w is a degree-4 "world" polynomial in rho / R_REF and m = Z / sqrt(X^2 + Y^2). The 2x2 affine stretch [[c, d], [e, 1]] is then applied to model sensor non-squareness / skew before translation by the principal point. This general polynomial radial profile (rather than a fixed functional form) is what lets OCam fit distortion curves other closed-form models cannot.

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates (u, v), origin top-left. On-axis points (sqrt(X^2+Y^2) < 1e-12) map to the principal point.
valid (ndarray, shape(...))

True iff the Newton solve for rho converged (|w(rho) + m*rho| < 1e-6) or the point is on-axis. See ds_msp.models.ocam_math.ocam_project.
References

Scaramuzza, D., Martinelli, A., Siegwart, R. "A Toolbox for Easily Calibrating Omnidirectional Cameras." IROS 2006.

Examples:

>>> import numpy as np
>>> from ds_msp.models import OCamModel
>>> m = OCamModel.sample()
>>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
>>> np.round(uv, 2)
array([[320., 240.]])
Source code in ds_msp/models/ocam.py 
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
def project(self, P):
    """Project camera-frame points via the Scaramuzza world polynomial.

    OCam has no closed-form projection: instead of a distortion formula
    applied to a normalized ray, it solves for the radial image distance
    ``rho`` satisfying ``w(rho) + m*rho = 0`` by Newton iteration, where
    ``w`` is a degree-4 "world" polynomial in ``rho / R_REF`` and
    ``m = Z / sqrt(X^2 + Y^2)``. The 2x2 affine stretch ``[[c, d], [e,
    1]]`` is then applied to model sensor non-squareness / skew before
    translation by the principal point. This general polynomial radial
    profile (rather than a fixed functional form) is what lets OCam fit
    distortion curves other closed-form models cannot.

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates ``(u, v)``, origin top-left.
        On-axis points (``sqrt(X^2+Y^2) < 1e-12``) map to the principal
        point.
    valid : ndarray, shape (...,)
        ``True`` iff the Newton solve for ``rho`` converged
        (``|w(rho) + m*rho| < 1e-6``) or the point is on-axis. See
        ``ds_msp.models.ocam_math.ocam_project``.

    References
    ----------
    Scaramuzza, D., Martinelli, A., Siegwart, R. "A Toolbox for Easily
    Calibrating Omnidirectional Cameras." IROS 2006.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import OCamModel
    >>> m = OCamModel.sample()
    >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
    >>> np.round(uv, 2)
    array([[320., 240.]])
    """
    u, v, valid = ocam_project(np.asarray(P, dtype=np.float64), *self._p())
    return np.stack([u, v], axis=-1), valid

project_jacobian 

project_jacobian(P)

Project with analytic derivatives via the implicit function theorem.

rho is defined implicitly by w(rho) + m*rho = 0; its derivatives w.r.t. both the 3D point and the polynomial coefficients are obtained without differentiating through the Newton loop, via the implicit function theorem (drho = -dF / (dF/drho)).

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates, identical to project.
J_point (ndarray, shape(..., 2, 3))

d(u, v) / d(X, Y, Z).
J_param (ndarray, shape(..., 2, 10))

d(u, v) / d(cx, cy, c, d, e, a0, a1, a2, a3, a4), columns in param_names order.
valid (ndarray, shape(...))

Projectability mask, identical condition to project.
References

Scaramuzza, D., Martinelli, A., Siegwart, R. IROS 2006 (Jacobian is this repo's own derivation via the implicit function theorem; verified by finite-difference check, relative error <= 1e-6, see pytest -m jac).

Examples:

>>> import numpy as np
>>> from ds_msp.models import OCamModel
>>> m = OCamModel.sample()
>>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
>>> J_point.shape, J_param.shape
((1, 2, 3), (1, 2, 10))
Source code in ds_msp/models/ocam.py 
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
def project_jacobian(self, P):
    """Project with analytic derivatives via the implicit function theorem.

    ``rho`` is defined implicitly by ``w(rho) + m*rho = 0``; its
    derivatives w.r.t. both the 3D point and the polynomial coefficients
    are obtained without differentiating through the Newton loop, via
    the implicit function theorem (``drho = -dF / (dF/drho)``).

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates, identical to `project`.
    J_point : ndarray, shape (..., 2, 3)
        ``d(u, v) / d(X, Y, Z)``.
    J_param : ndarray, shape (..., 2, 10)
        ``d(u, v) / d(cx, cy, c, d, e, a0, a1, a2, a3, a4)``, columns in
        `param_names` order.
    valid : ndarray, shape (...,)
        Projectability mask, identical condition to `project`.

    References
    ----------
    Scaramuzza, D., Martinelli, A., Siegwart, R. IROS 2006 (Jacobian is
    this repo's own derivation via the implicit function theorem;
    verified by finite-difference check, relative error <= 1e-6, see
    ``pytest -m jac``).

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import OCamModel
    >>> m = OCamModel.sample()
    >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
    >>> J_point.shape, J_param.shape
    ((1, 2, 3), (1, 2, 10))
    """
    u, v, J_point, J_param, valid = ocam_project_jacobian(
        np.asarray(P, dtype=np.float64), *self._p())
    return np.stack([u, v], axis=-1), J_point, J_param, valid

sample classmethod 

sample() -> 'OCamModel'

Realistic instance for contract testing (a synthetic omni-camera).

Polynomial argument is normalized by R_REF=100, so a2..a4 are O(1).

Source code in ds_msp/models/ocam.py 
29
30
31
32
33
34
35
36
@classmethod
def sample(cls) -> "OCamModel":
    """Realistic instance for contract testing (a synthetic omni-camera).

    Polynomial argument is normalized by ``R_REF=100``, so ``a2..a4``
    are O(1).
    """
    return cls(320.0, 240.0, 1.0, 0.0, 0.0, -220.0, 0.0, 6.0, -1.5, 0.25)

to_dict 

to_dict() -> dict

Serialize to {"model": "ocam", "cx": ..., ..., "a4": ...}.

Source code in ds_msp/models/ocam.py 
214
215
216
217
218
def to_dict(self) -> dict:
    """Serialize to ``{"model": "ocam", "cx": ..., ..., "a4": ...}``."""
    d = {"model": self.name}
    d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
    return d

unproject 

unproject(uv)

Unproject pixels to unit bearing rays (direct polynomial evaluation).

Inverts the affine stretch linearly, then evaluates the world polynomial w(rho) directly (no root-finding needed in this direction) to form ray = normalize([x, y, -w(rho)]). See ds_msp.models.ocam_math.ocam_unproject.

Parameters:

Name Type Description Default
uv (ndarray, shape(..., 2))

Pixel coordinates (u, v), origin top-left.

 required

Returns:

Name Type Description
rays (ndarray, shape(..., 3))

Unit-norm camera-frame bearing vectors, +Z forward.
valid (ndarray, shape(...))

True iff the resulting ray is finite. Invalid rays are zeroed, never NaN.

Examples:

>>> import numpy as np
>>> from ds_msp.models import OCamModel
>>> m = OCamModel.sample()
>>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
>>> np.round(rays, 4)
array([[0., 0., 1.]])
Source code in ds_msp/models/ocam.py 
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
def unproject(self, uv):
    """Unproject pixels to unit bearing rays (direct polynomial evaluation).

    Inverts the affine stretch linearly, then evaluates the world
    polynomial ``w(rho)`` directly (no root-finding needed in this
    direction) to form ``ray = normalize([x, y, -w(rho)])``. See
    ``ds_msp.models.ocam_math.ocam_unproject``.

    Parameters
    ----------
    uv : ndarray, shape (..., 2)
        Pixel coordinates ``(u, v)``, origin top-left.

    Returns
    -------
    rays : ndarray, shape (..., 3)
        Unit-norm camera-frame bearing vectors, +Z forward.
    valid : ndarray, shape (...,)
        ``True`` iff the resulting ray is finite. Invalid rays are
        zeroed, never NaN.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import OCamModel
    >>> m = OCamModel.sample()
    >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
    >>> np.round(rays, 4)
    array([[0., 0., 1.]])
    """
    return ocam_unproject(np.asarray(uv, dtype=np.float64), *self._p())

Radial-Tangential (RadTan / OpenCV pinhole)

ds_msp.models.RadTanModel

Pinhole with radial-tangential distortion. Satisfies CameraModel.

distortion returns [k1, k2, p1, p2, k3] (OpenCV order) for direct use with cv2.projectPoints / cv2.undistortPoints. Narrow-FOV: only models rays with z > 0.

Source code in ds_msp/models/radtan.py 
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
class RadTanModel:
    """Pinhole with radial-tangential distortion. Satisfies ``CameraModel``.

    ``distortion`` returns ``[k1, k2, p1, p2, k3]`` (OpenCV order) for direct use
    with ``cv2.projectPoints`` / ``cv2.undistortPoints``. Narrow-FOV: only models
    rays with ``z > 0``.
    """

    name: ClassVar[str] = "radtan"
    param_names: ClassVar[Tuple[str, ...]] = (
        "fx", "fy", "cx", "cy", "k1", "k2", "p1", "p2", "k3")

    def __init__(self, fx, fy, cx, cy, k1=0.0, k2=0.0, p1=0.0, p2=0.0, k3=0.0) -> None:
        self.fx, self.fy, self.cx, self.cy = float(fx), float(fy), float(cx), float(cy)
        self.k1, self.k2, self.k3 = float(k1), float(k2), float(k3)
        self.p1, self.p2 = float(p1), float(p2)

    @classmethod
    def sample(cls) -> "RadTanModel":
        """Realistic instance for contract testing (a narrow-FOV pinhole lens)."""
        return cls(600.0, 602.0, 320.0, 240.0, -0.12, 0.05, 0.001, -0.0015, 0.008)

    @property
    def params(self) -> np.ndarray:
        """Flat parameter vector ``[fx, fy, cx, cy, k1, k2, p1, p2, k3]``."""
        return np.array([self.fx, self.fy, self.cx, self.cy,
                         self.k1, self.k2, self.p1, self.p2, self.k3], dtype=np.float64)

    @property
    def K(self) -> np.ndarray:
        """3x3 pinhole intrinsic matrix built from ``fx, fy, cx, cy``."""
        return np.array([[self.fx, 0.0, self.cx], [0.0, self.fy, self.cy],
                         [0.0, 0.0, 1.0]], dtype=np.float64)

    @property
    def distortion(self) -> np.ndarray:
        """OpenCV distCoeffs order [k1, k2, p1, p2, k3]."""
        return np.array([self.k1, self.k2, self.p1, self.p2, self.k3], dtype=np.float64)

    def _coeffs(self):
        return (self.k1, self.k2, self.p1, self.p2, self.k3)

    def project(self, P):
        """Project camera-frame points via the classic Brown-Conrady model.

        The pinhole projects onto the normalized plane (``a=x/z, b=y/z``)
        and applies polynomial radial distortion (``k1, k2, k3`` on
        ``r2=a^2+b^2``) plus tangential distortion (``p1, p2``, modeling lens
        decentering). Unlike the other 7 models in this package, RadTan does
        **not** attempt to represent fisheye FOV: it is only valid for
        ``z > 0`` and cannot represent rays at or beyond 90° off-axis.

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates ``(u, v)``, origin top-left.
        valid : ndarray, shape (...,)
            ``True`` iff ``z > 1e-9``. Points with ``z <= 0`` are masked and
            their pixels zeroed. See ``ds_msp.models.radtan_math.radtan_project``.

        References
        ----------
        Brown, D. C. "Decentering Distortion of Lenses." Photogrammetric
        Engineering, 1966 (tangential terms); OpenCV
        ``distCoeffs = [k1, k2, p1, p2, k3]`` convention, see
        ``cv2.projectPoints``.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import RadTanModel
        >>> m = RadTanModel.sample()
        >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
        >>> np.round(uv, 2)
        array([[320., 240.]])
        """
        u, v, valid = radtan_project(np.asarray(P, dtype=np.float64),
                                     self.fx, self.fy, self.cx, self.cy, *self._coeffs())
        return np.stack([u, v], axis=-1), valid

    def unproject(self, uv):
        """Unproject pixels to unit bearing rays via fixed-point distortion inversion.

        Runs 20 fixed-point iterations of the standard OpenCV-style inverse
        (no closed form exists for the Brown-Conrady distortion), then
        re-distorts the recovered normalized coordinates to check
        consistency with the input.

        Parameters
        ----------
        uv : ndarray, shape (..., 2)
            Pixel coordinates ``(u, v)``, origin top-left.

        Returns
        -------
        rays : ndarray, shape (..., 3)
            Unit-norm camera-frame bearing vectors, +Z forward.
        valid : ndarray, shape (...,)
            ``True`` iff re-distorting the recovered point reproduces the
            input to within ``1e-6``. Invalid rays are zeroed, never NaN.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import RadTanModel
        >>> m = RadTanModel.sample()
        >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
        >>> np.round(rays, 4)
        array([[0., 0., 1.]])
        """
        return radtan_unproject(np.asarray(uv, dtype=np.float64),
                                self.fx, self.fy, self.cx, self.cy, *self._coeffs())

    def project_jacobian(self, P):
        """Project with analytic derivatives (no autodiff, no finite differences).

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates, identical to `project`.
        J_point : ndarray, shape (..., 2, 3)
            ``d(u, v) / d(x, y, z)``.
        J_param : ndarray, shape (..., 2, 9)
            ``d(u, v) / d(fx, fy, cx, cy, k1, k2, p1, p2, k3)``, columns in
            `param_names` order.
        valid : ndarray, shape (...,)
            Projectability mask, identical condition to `project`.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import RadTanModel
        >>> m = RadTanModel.sample()
        >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
        >>> J_point.shape, J_param.shape
        ((1, 2, 3), (1, 2, 9))
        """
        u, v, J_point, J_param, valid = radtan_project_jacobian(
            np.asarray(P, dtype=np.float64),
            self.fx, self.fy, self.cx, self.cy, *self._coeffs())
        return np.stack([u, v], axis=-1), J_point, J_param, valid

    @classmethod
    def from_params(cls, p):
        """Build from a flat ``[fx, fy, cx, cy, k1, k2, p1, p2, k3]`` vector."""
        return cls(*np.asarray(p, dtype=np.float64).ravel())

    @classmethod
    def param_bounds(cls):
        """Optimizer bounds: radial/tangential terms modestly bounded, focal/center wide-open."""
        lb = np.array([1.0, 1.0, -1e5, -1e5, -2.0, -2.0, -1.0, -1.0, -2.0], dtype=np.float64)
        ub = np.array([1e5, 1e5, 1e5, 1e5, 2.0, 2.0, 1.0, 1.0, 2.0], dtype=np.float64)
        return lb, ub

    def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
        """Seed ``fx,fy,cx,cy`` from `K_seed`; zero all distortion coefficients."""
        self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
        self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
        self.k1 = self.k2 = self.k3 = self.p1 = self.p2 = 0.0

    def to_dict(self) -> dict:
        """Serialize to ``{"model": "radtan", "fx": ..., ..., "k3": ...}``."""
        d = {"model": self.name}
        d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
        return d

    @classmethod
    def from_dict(cls, d: dict) -> "RadTanModel":
        """Reconstruct from :meth:`to_dict` output."""
        return cls(**{k: d[k] for k in cls.param_names})

    def __repr__(self) -> str:
        return ("RadTanModel(fx={:.3f}, fy={:.3f}, cx={:.3f}, cy={:.3f}, "
                "k1={:.5f}, k2={:.5f}, p1={:.5f}, p2={:.5f}, k3={:.5f})").format(
                    self.fx, self.fy, self.cx, self.cy,
                    self.k1, self.k2, self.p1, self.p2, self.k3)

K property 

K: ndarray

3x3 pinhole intrinsic matrix built from fx, fy, cx, cy.

distortion property 

distortion: ndarray

OpenCV distCoeffs order [k1, k2, p1, p2, k3].

params property 

params: ndarray

Flat parameter vector [fx, fy, cx, cy, k1, k2, p1, p2, k3].

from_dict classmethod 

from_dict(d: dict) -> 'RadTanModel'

Reconstruct from :meth:to_dict output.

Source code in ds_msp/models/radtan.py 
188
189
190
191
@classmethod
def from_dict(cls, d: dict) -> "RadTanModel":
    """Reconstruct from :meth:`to_dict` output."""
    return cls(**{k: d[k] for k in cls.param_names})

from_params classmethod 

from_params(p)

Build from a flat [fx, fy, cx, cy, k1, k2, p1, p2, k3] vector.

Source code in ds_msp/models/radtan.py 
164
165
166
167
@classmethod
def from_params(cls, p):
    """Build from a flat ``[fx, fy, cx, cy, k1, k2, p1, p2, k3]`` vector."""
    return cls(*np.asarray(p, dtype=np.float64).ravel())

initialize_from_correspondences 

initialize_from_correspondences(K_seed, rays, pixels) -> None

Seed fx,fy,cx,cy from K_seed; zero all distortion coefficients.

Source code in ds_msp/models/radtan.py 
176
177
178
179
180
def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
    """Seed ``fx,fy,cx,cy`` from `K_seed`; zero all distortion coefficients."""
    self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
    self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
    self.k1 = self.k2 = self.k3 = self.p1 = self.p2 = 0.0

param_bounds classmethod 

param_bounds()

Optimizer bounds: radial/tangential terms modestly bounded, focal/center wide-open.

Source code in ds_msp/models/radtan.py 
169
170
171
172
173
174
@classmethod
def param_bounds(cls):
    """Optimizer bounds: radial/tangential terms modestly bounded, focal/center wide-open."""
    lb = np.array([1.0, 1.0, -1e5, -1e5, -2.0, -2.0, -1.0, -1.0, -2.0], dtype=np.float64)
    ub = np.array([1e5, 1e5, 1e5, 1e5, 2.0, 2.0, 1.0, 1.0, 2.0], dtype=np.float64)
    return lb, ub

project 

project(P)

Project camera-frame points via the classic Brown-Conrady model.

The pinhole projects onto the normalized plane (a=x/z, b=y/z) and applies polynomial radial distortion (k1, k2, k3 on r2=a^2+b^2) plus tangential distortion (p1, p2, modeling lens decentering). Unlike the other 7 models in this package, RadTan does not attempt to represent fisheye FOV: it is only valid for z > 0 and cannot represent rays at or beyond 90° off-axis.

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates (u, v), origin top-left.
valid (ndarray, shape(...))

True iff z > 1e-9. Points with z <= 0 are masked and their pixels zeroed. See ds_msp.models.radtan_math.radtan_project.
References

Brown, D. C. "Decentering Distortion of Lenses." Photogrammetric Engineering, 1966 (tangential terms); OpenCV distCoeffs = [k1, k2, p1, p2, k3] convention, see cv2.projectPoints.

Examples:

>>> import numpy as np
>>> from ds_msp.models import RadTanModel
>>> m = RadTanModel.sample()
>>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
>>> np.round(uv, 2)
array([[320., 240.]])
Source code in ds_msp/models/radtan.py 
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
def project(self, P):
    """Project camera-frame points via the classic Brown-Conrady model.

    The pinhole projects onto the normalized plane (``a=x/z, b=y/z``)
    and applies polynomial radial distortion (``k1, k2, k3`` on
    ``r2=a^2+b^2``) plus tangential distortion (``p1, p2``, modeling lens
    decentering). Unlike the other 7 models in this package, RadTan does
    **not** attempt to represent fisheye FOV: it is only valid for
    ``z > 0`` and cannot represent rays at or beyond 90° off-axis.

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates ``(u, v)``, origin top-left.
    valid : ndarray, shape (...,)
        ``True`` iff ``z > 1e-9``. Points with ``z <= 0`` are masked and
        their pixels zeroed. See ``ds_msp.models.radtan_math.radtan_project``.

    References
    ----------
    Brown, D. C. "Decentering Distortion of Lenses." Photogrammetric
    Engineering, 1966 (tangential terms); OpenCV
    ``distCoeffs = [k1, k2, p1, p2, k3]`` convention, see
    ``cv2.projectPoints``.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import RadTanModel
    >>> m = RadTanModel.sample()
    >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
    >>> np.round(uv, 2)
    array([[320., 240.]])
    """
    u, v, valid = radtan_project(np.asarray(P, dtype=np.float64),
                                 self.fx, self.fy, self.cx, self.cy, *self._coeffs())
    return np.stack([u, v], axis=-1), valid

project_jacobian 

project_jacobian(P)

Project with analytic derivatives (no autodiff, no finite differences).

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates, identical to project.
J_point (ndarray, shape(..., 2, 3))

d(u, v) / d(x, y, z).
J_param (ndarray, shape(..., 2, 9))

d(u, v) / d(fx, fy, cx, cy, k1, k2, p1, p2, k3), columns in param_names order.
valid (ndarray, shape(...))

Projectability mask, identical condition to project.

Examples:

>>> import numpy as np
>>> from ds_msp.models import RadTanModel
>>> m = RadTanModel.sample()
>>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
>>> J_point.shape, J_param.shape
((1, 2, 3), (1, 2, 9))
Source code in ds_msp/models/radtan.py 
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
def project_jacobian(self, P):
    """Project with analytic derivatives (no autodiff, no finite differences).

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates, identical to `project`.
    J_point : ndarray, shape (..., 2, 3)
        ``d(u, v) / d(x, y, z)``.
    J_param : ndarray, shape (..., 2, 9)
        ``d(u, v) / d(fx, fy, cx, cy, k1, k2, p1, p2, k3)``, columns in
        `param_names` order.
    valid : ndarray, shape (...,)
        Projectability mask, identical condition to `project`.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import RadTanModel
    >>> m = RadTanModel.sample()
    >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
    >>> J_point.shape, J_param.shape
    ((1, 2, 3), (1, 2, 9))
    """
    u, v, J_point, J_param, valid = radtan_project_jacobian(
        np.asarray(P, dtype=np.float64),
        self.fx, self.fy, self.cx, self.cy, *self._coeffs())
    return np.stack([u, v], axis=-1), J_point, J_param, valid

sample classmethod 

sample() -> 'RadTanModel'

Realistic instance for contract testing (a narrow-FOV pinhole lens).

Source code in ds_msp/models/radtan.py 
29
30
31
32
@classmethod
def sample(cls) -> "RadTanModel":
    """Realistic instance for contract testing (a narrow-FOV pinhole lens)."""
    return cls(600.0, 602.0, 320.0, 240.0, -0.12, 0.05, 0.001, -0.0015, 0.008)

to_dict 

to_dict() -> dict

Serialize to {"model": "radtan", "fx": ..., ..., "k3": ...}.

Source code in ds_msp/models/radtan.py 
182
183
184
185
186
def to_dict(self) -> dict:
    """Serialize to ``{"model": "radtan", "fx": ..., ..., "k3": ...}``."""
    d = {"model": self.name}
    d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
    return d

unproject 

unproject(uv)

Unproject pixels to unit bearing rays via fixed-point distortion inversion.

Runs 20 fixed-point iterations of the standard OpenCV-style inverse (no closed form exists for the Brown-Conrady distortion), then re-distorts the recovered normalized coordinates to check consistency with the input.

Parameters:

Name Type Description Default
uv (ndarray, shape(..., 2))

Pixel coordinates (u, v), origin top-left.

 required

Returns:

Name Type Description
rays (ndarray, shape(..., 3))

Unit-norm camera-frame bearing vectors, +Z forward.
valid (ndarray, shape(...))

True iff re-distorting the recovered point reproduces the input to within 1e-6. Invalid rays are zeroed, never NaN.

Examples:

>>> import numpy as np
>>> from ds_msp.models import RadTanModel
>>> m = RadTanModel.sample()
>>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
>>> np.round(rays, 4)
array([[0., 0., 1.]])
Source code in ds_msp/models/radtan.py 
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
def unproject(self, uv):
    """Unproject pixels to unit bearing rays via fixed-point distortion inversion.

    Runs 20 fixed-point iterations of the standard OpenCV-style inverse
    (no closed form exists for the Brown-Conrady distortion), then
    re-distorts the recovered normalized coordinates to check
    consistency with the input.

    Parameters
    ----------
    uv : ndarray, shape (..., 2)
        Pixel coordinates ``(u, v)``, origin top-left.

    Returns
    -------
    rays : ndarray, shape (..., 3)
        Unit-norm camera-frame bearing vectors, +Z forward.
    valid : ndarray, shape (...,)
        ``True`` iff re-distorting the recovered point reproduces the
        input to within ``1e-6``. Invalid rays are zeroed, never NaN.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import RadTanModel
    >>> m = RadTanModel.sample()
    >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
    >>> np.round(rays, 4)
    array([[0., 0., 1.]])
    """
    return radtan_unproject(np.asarray(uv, dtype=np.float64),
                            self.fx, self.fy, self.cx, self.cy, *self._coeffs())

Unified Camera Model (UCM)

ds_msp.models.UCMModel

Unified Camera Model (single sphere). Satisfies CameraModel.

Source code in ds_msp/models/ucm.py 
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
class UCMModel:
    """Unified Camera Model (single sphere). Satisfies ``CameraModel``."""

    name: ClassVar[str] = "ucm"
    param_names: ClassVar[Tuple[str, ...]] = ("fx", "fy", "cx", "cy", "alpha")

    def __init__(self, fx: float, fy: float, cx: float, cy: float, alpha: float) -> None:
        self.fx = float(fx)
        self.fy = float(fy)
        self.cx = float(cx)
        self.cy = float(cy)
        self.alpha = float(alpha)

    @classmethod
    def sample(cls) -> "UCMModel":
        """Realistic instance for contract testing (the bundled calibration)."""
        return cls(711.57, 711.24, 949.18, 518.81, 0.62)

    @property
    def params(self) -> np.ndarray:
        """Flat parameter vector ``[fx, fy, cx, cy, alpha]``."""
        return np.array([self.fx, self.fy, self.cx, self.cy, self.alpha], dtype=np.float64)

    @property
    def K(self) -> np.ndarray:
        """3x3 pinhole intrinsic matrix built from ``fx, fy, cx, cy``."""
        return np.array([[self.fx, 0.0, self.cx], [0.0, self.fy, self.cy],
                         [0.0, 0.0, 1.0]], dtype=np.float64)

    @property
    def distortion(self) -> np.ndarray:
        """Distortion tail ``[alpha]`` (perspective blend between sphere and plane)."""
        return np.array([self.alpha], dtype=np.float64)

    def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Project camera-frame points via the single unit-sphere composition.

        UCM re-centers the point onto a single unit sphere (radius
        ``d = sqrt(x^2+y^2+z^2)``), then perspective-divides from a point
        blended between the sphere's surface and the pinhole plane by
        ``alpha``: ``den = alpha*d + (1-alpha)*z``. It is the one-sphere,
        one-parameter special case of `DoubleSphereModel` (``xi=0``).

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates ``(u, v)``, origin top-left.
        valid : ndarray, shape (...,)
            ``True`` iff the point lies in the tilted half-space
            ``z > -w(alpha) * d`` (*not* the naive ``z > 0``) and the
            perspective-division denominator is bounded away from zero. See
            ``ds_msp.models.ucm_math.ucm_project``.

        References
        ----------
        Geyer, C., Daniilidis, K. "A Unifying Theory for Central Panoramic
        Systems." ECCV 2000; Mei, C., Rives, P. "Single View Point
        Omnidirectional Camera Calibration from Planar Grids." ICRA 2007.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import UCMModel
        >>> m = UCMModel.sample()
        >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
        >>> np.round(uv, 2)
        array([[949.18, 518.81]])
        """
        u, v, valid = ucm_project(np.asarray(P, dtype=np.float64),
                                  self.fx, self.fy, self.cx, self.cy, self.alpha)
        return np.stack([u, v], axis=-1), valid

    def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
        """Unproject pixels to unit bearing rays (closed form).

        Parameters
        ----------
        uv : ndarray, shape (..., 2)
            Pixel coordinates ``(u, v)``, origin top-left.

        Returns
        -------
        rays : ndarray, shape (..., 3)
            Unit-norm camera-frame bearing vectors, +Z forward.
        valid : ndarray, shape (...,)
            ``True`` iff the pixel lies within the sphere's reachable disc
            and the perspective-division denominator is bounded away from
            zero. Invalid rays are zeroed, never NaN. See
            ``ds_msp.models.ucm_math.ucm_unproject``.

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import UCMModel
        >>> m = UCMModel.sample()
        >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
        >>> np.round(rays, 4)
        array([[0., 0., 1.]])
        """
        return ucm_unproject(np.asarray(uv, dtype=np.float64),
                             self.fx, self.fy, self.cx, self.cy, self.alpha)

    def project_jacobian(self, P):
        """Project with analytic derivatives (no autodiff, no finite differences).

        Parameters
        ----------
        P : ndarray, shape (..., 3)
            Camera-frame points (meters), +Z forward.

        Returns
        -------
        uv : ndarray, shape (..., 2)
            Projected pixel coordinates, identical to `project`.
        J_point : ndarray, shape (..., 2, 3)
            ``d(u, v) / d(x, y, z)``.
        J_param : ndarray, shape (..., 2, 5)
            ``d(u, v) / d(fx, fy, cx, cy, alpha)``, columns in `param_names`
            order.
        valid : ndarray, shape (...,)
            Projectability mask, identical condition to `project`.

        References
        ----------
        Geyer, C., Daniilidis, K. ECCV 2000; Mei, C., Rives, P. ICRA 2007
        (closed-form Jacobian derived from the forward map; verified here by
        finite-difference check, relative error <= 1e-6, see ``pytest -m jac``).

        Examples
        --------
        >>> import numpy as np
        >>> from ds_msp.models import UCMModel
        >>> m = UCMModel.sample()
        >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
        >>> J_point.shape, J_param.shape
        ((1, 2, 3), (1, 2, 5))
        """
        u, v, J_point, J_param, valid = ucm_project_jacobian(
            np.asarray(P, dtype=np.float64),
            self.fx, self.fy, self.cx, self.cy, self.alpha)
        return np.stack([u, v], axis=-1), J_point, J_param, valid

    @classmethod
    def from_params(cls, p: np.ndarray) -> "UCMModel":
        """Build from a flat ``[fx, fy, cx, cy, alpha]`` vector."""
        return cls(*np.asarray(p, dtype=np.float64).ravel())

    @classmethod
    def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
        """Optimizer bounds: ``alpha in (0, 1)``, focal/center wide-open."""
        lb = np.array([1.0, 1.0, -1e5, -1e5, 1e-6], dtype=np.float64)
        ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0 - 1e-6], dtype=np.float64)
        return lb, ub

    def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
        """Seed ``fx,fy,cx,cy`` from `K_seed`; solve ``alpha`` by linear least squares
        from unit-ray/pixel correspondences."""
        self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
        self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
        rays = np.asarray(rays, dtype=np.float64)
        x, y, z = rays[:, 0], rays[:, 1], rays[:, 2]
        mx = (pixels[:, 0] - self.cx) / self.fx
        my = (pixels[:, 1] - self.cy) / self.fy
        # unit rays (d = 1): alpha = (x - mx*z) / (mx*(1 - z)), solved linearly.
        A = np.concatenate([mx * (1.0 - z), my * (1.0 - z)])
        b = np.concatenate([x - mx * z, y - my * z])
        denom = float(A @ A)
        self.alpha = float(np.clip((A @ b) / denom, 1e-6, 1.0 - 1e-6)) if denom > 1e-12 else 0.5

    def to_dict(self) -> dict:
        """Serialize to ``{"model": "ucm", "fx": ..., ..., "alpha": ...}``."""
        d = {"model": self.name}
        d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
        return d

    @classmethod
    def from_dict(cls, d: dict) -> "UCMModel":
        """Reconstruct from :meth:`to_dict` output."""
        return cls(**{k: d[k] for k in cls.param_names})

    def __repr__(self) -> str:
        return "UCMModel(fx={:.3f}, fy={:.3f}, cx={:.3f}, cy={:.3f}, alpha={:.4f})".format(
            self.fx, self.fy, self.cx, self.cy, self.alpha)

K property 

K: ndarray

3x3 pinhole intrinsic matrix built from fx, fy, cx, cy.

distortion property 

distortion: ndarray

Distortion tail [alpha] (perspective blend between sphere and plane).

params property 

params: ndarray

Flat parameter vector [fx, fy, cx, cy, alpha].

from_dict classmethod 

from_dict(d: dict) -> 'UCMModel'

Reconstruct from :meth:to_dict output.

Source code in ds_msp/models/ucm.py 
192
193
194
195
@classmethod
def from_dict(cls, d: dict) -> "UCMModel":
    """Reconstruct from :meth:`to_dict` output."""
    return cls(**{k: d[k] for k in cls.param_names})

from_params classmethod 

from_params(p: ndarray) -> 'UCMModel'

Build from a flat [fx, fy, cx, cy, alpha] vector.

Source code in ds_msp/models/ucm.py 
159
160
161
162
@classmethod
def from_params(cls, p: np.ndarray) -> "UCMModel":
    """Build from a flat ``[fx, fy, cx, cy, alpha]`` vector."""
    return cls(*np.asarray(p, dtype=np.float64).ravel())

initialize_from_correspondences 

initialize_from_correspondences(K_seed, rays, pixels) -> None

Seed fx,fy,cx,cy from K_seed; solve alpha by linear least squares from unit-ray/pixel correspondences.

Source code in ds_msp/models/ucm.py 
171
172
173
174
175
176
177
178
179
180
181
182
183
184
def initialize_from_correspondences(self, K_seed, rays, pixels) -> None:
    """Seed ``fx,fy,cx,cy`` from `K_seed`; solve ``alpha`` by linear least squares
    from unit-ray/pixel correspondences."""
    self.fx, self.fy = float(K_seed[0, 0]), float(K_seed[1, 1])
    self.cx, self.cy = float(K_seed[0, 2]), float(K_seed[1, 2])
    rays = np.asarray(rays, dtype=np.float64)
    x, y, z = rays[:, 0], rays[:, 1], rays[:, 2]
    mx = (pixels[:, 0] - self.cx) / self.fx
    my = (pixels[:, 1] - self.cy) / self.fy
    # unit rays (d = 1): alpha = (x - mx*z) / (mx*(1 - z)), solved linearly.
    A = np.concatenate([mx * (1.0 - z), my * (1.0 - z)])
    b = np.concatenate([x - mx * z, y - my * z])
    denom = float(A @ A)
    self.alpha = float(np.clip((A @ b) / denom, 1e-6, 1.0 - 1e-6)) if denom > 1e-12 else 0.5

param_bounds classmethod 

param_bounds() -> Tuple[np.ndarray, np.ndarray]

Optimizer bounds: alpha in (0, 1), focal/center wide-open.

Source code in ds_msp/models/ucm.py 
164
165
166
167
168
169
@classmethod
def param_bounds(cls) -> Tuple[np.ndarray, np.ndarray]:
    """Optimizer bounds: ``alpha in (0, 1)``, focal/center wide-open."""
    lb = np.array([1.0, 1.0, -1e5, -1e5, 1e-6], dtype=np.float64)
    ub = np.array([1e5, 1e5, 1e5, 1e5, 1.0 - 1e-6], dtype=np.float64)
    return lb, ub

project 

project(P: ndarray) -> Tuple[np.ndarray, np.ndarray]

Project camera-frame points via the single unit-sphere composition.

UCM re-centers the point onto a single unit sphere (radius d = sqrt(x^2+y^2+z^2)), then perspective-divides from a point blended between the sphere's surface and the pinhole plane by alpha: den = alpha*d + (1-alpha)*z. It is the one-sphere, one-parameter special case of DoubleSphereModel (xi=0).

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates (u, v), origin top-left.
valid (ndarray, shape(...))

True iff the point lies in the tilted half-space z > -w(alpha) * d (not the naive z > 0) and the perspective-division denominator is bounded away from zero. See ds_msp.models.ucm_math.ucm_project.
References

Geyer, C., Daniilidis, K. "A Unifying Theory for Central Panoramic Systems." ECCV 2000; Mei, C., Rives, P. "Single View Point Omnidirectional Camera Calibration from Planar Grids." ICRA 2007.

Examples:

>>> import numpy as np
>>> from ds_msp.models import UCMModel
>>> m = UCMModel.sample()
>>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
>>> np.round(uv, 2)
array([[949.18, 518.81]])
Source code in ds_msp/models/ucm.py 
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
def project(self, P: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Project camera-frame points via the single unit-sphere composition.

    UCM re-centers the point onto a single unit sphere (radius
    ``d = sqrt(x^2+y^2+z^2)``), then perspective-divides from a point
    blended between the sphere's surface and the pinhole plane by
    ``alpha``: ``den = alpha*d + (1-alpha)*z``. It is the one-sphere,
    one-parameter special case of `DoubleSphereModel` (``xi=0``).

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates ``(u, v)``, origin top-left.
    valid : ndarray, shape (...,)
        ``True`` iff the point lies in the tilted half-space
        ``z > -w(alpha) * d`` (*not* the naive ``z > 0``) and the
        perspective-division denominator is bounded away from zero. See
        ``ds_msp.models.ucm_math.ucm_project``.

    References
    ----------
    Geyer, C., Daniilidis, K. "A Unifying Theory for Central Panoramic
    Systems." ECCV 2000; Mei, C., Rives, P. "Single View Point
    Omnidirectional Camera Calibration from Planar Grids." ICRA 2007.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import UCMModel
    >>> m = UCMModel.sample()
    >>> uv, valid = m.project(np.array([[0.0, 0.0, 1.0]]))
    >>> np.round(uv, 2)
    array([[949.18, 518.81]])
    """
    u, v, valid = ucm_project(np.asarray(P, dtype=np.float64),
                              self.fx, self.fy, self.cx, self.cy, self.alpha)
    return np.stack([u, v], axis=-1), valid

project_jacobian 

project_jacobian(P)

Project with analytic derivatives (no autodiff, no finite differences).

Parameters:

Name Type Description Default
P (ndarray, shape(..., 3))

Camera-frame points (meters), +Z forward.

 required

Returns:

Name Type Description
uv (ndarray, shape(..., 2))

Projected pixel coordinates, identical to project.
J_point (ndarray, shape(..., 2, 3))

d(u, v) / d(x, y, z).
J_param (ndarray, shape(..., 2, 5))

d(u, v) / d(fx, fy, cx, cy, alpha), columns in param_names order.
valid (ndarray, shape(...))

Projectability mask, identical condition to project.
References

Geyer, C., Daniilidis, K. ECCV 2000; Mei, C., Rives, P. ICRA 2007 (closed-form Jacobian derived from the forward map; verified here by finite-difference check, relative error <= 1e-6, see pytest -m jac).

Examples:

>>> import numpy as np
>>> from ds_msp.models import UCMModel
>>> m = UCMModel.sample()
>>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
>>> J_point.shape, J_param.shape
((1, 2, 3), (1, 2, 5))
Source code in ds_msp/models/ucm.py 
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
def project_jacobian(self, P):
    """Project with analytic derivatives (no autodiff, no finite differences).

    Parameters
    ----------
    P : ndarray, shape (..., 3)
        Camera-frame points (meters), +Z forward.

    Returns
    -------
    uv : ndarray, shape (..., 2)
        Projected pixel coordinates, identical to `project`.
    J_point : ndarray, shape (..., 2, 3)
        ``d(u, v) / d(x, y, z)``.
    J_param : ndarray, shape (..., 2, 5)
        ``d(u, v) / d(fx, fy, cx, cy, alpha)``, columns in `param_names`
        order.
    valid : ndarray, shape (...,)
        Projectability mask, identical condition to `project`.

    References
    ----------
    Geyer, C., Daniilidis, K. ECCV 2000; Mei, C., Rives, P. ICRA 2007
    (closed-form Jacobian derived from the forward map; verified here by
    finite-difference check, relative error <= 1e-6, see ``pytest -m jac``).

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import UCMModel
    >>> m = UCMModel.sample()
    >>> uv, J_point, J_param, valid = m.project_jacobian(np.array([[0.0, 0.0, 1.0]]))
    >>> J_point.shape, J_param.shape
    ((1, 2, 3), (1, 2, 5))
    """
    u, v, J_point, J_param, valid = ucm_project_jacobian(
        np.asarray(P, dtype=np.float64),
        self.fx, self.fy, self.cx, self.cy, self.alpha)
    return np.stack([u, v], axis=-1), J_point, J_param, valid

sample classmethod 

sample() -> 'UCMModel'

Realistic instance for contract testing (the bundled calibration).

Source code in ds_msp/models/ucm.py 
25
26
27
28
@classmethod
def sample(cls) -> "UCMModel":
    """Realistic instance for contract testing (the bundled calibration)."""
    return cls(711.57, 711.24, 949.18, 518.81, 0.62)

to_dict 

to_dict() -> dict

Serialize to {"model": "ucm", "fx": ..., ..., "alpha": ...}.

Source code in ds_msp/models/ucm.py 
186
187
188
189
190
def to_dict(self) -> dict:
    """Serialize to ``{"model": "ucm", "fx": ..., ..., "alpha": ...}``."""
    d = {"model": self.name}
    d.update({k: float(v) for k, v in zip(self.param_names, self.params)})
    return d

unproject 

unproject(uv: ndarray) -> Tuple[np.ndarray, np.ndarray]

Unproject pixels to unit bearing rays (closed form).

Parameters:

Name Type Description Default
uv (ndarray, shape(..., 2))

Pixel coordinates (u, v), origin top-left.

 required

Returns:

Name Type Description
rays (ndarray, shape(..., 3))

Unit-norm camera-frame bearing vectors, +Z forward.
valid (ndarray, shape(...))

True iff the pixel lies within the sphere's reachable disc and the perspective-division denominator is bounded away from zero. Invalid rays are zeroed, never NaN. See ds_msp.models.ucm_math.ucm_unproject.

Examples:

>>> import numpy as np
>>> from ds_msp.models import UCMModel
>>> m = UCMModel.sample()
>>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
>>> np.round(rays, 4)
array([[0., 0., 1.]])
Source code in ds_msp/models/ucm.py 
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
def unproject(self, uv: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Unproject pixels to unit bearing rays (closed form).

    Parameters
    ----------
    uv : ndarray, shape (..., 2)
        Pixel coordinates ``(u, v)``, origin top-left.

    Returns
    -------
    rays : ndarray, shape (..., 3)
        Unit-norm camera-frame bearing vectors, +Z forward.
    valid : ndarray, shape (...,)
        ``True`` iff the pixel lies within the sphere's reachable disc
        and the perspective-division denominator is bounded away from
        zero. Invalid rays are zeroed, never NaN. See
        ``ds_msp.models.ucm_math.ucm_unproject``.

    Examples
    --------
    >>> import numpy as np
    >>> from ds_msp.models import UCMModel
    >>> m = UCMModel.sample()
    >>> rays, valid = m.unproject(np.array([[m.cx, m.cy]]))
    >>> np.round(rays, 4)
    array([[0., 0., 1.]])
    """
    return ucm_unproject(np.asarray(uv, dtype=np.float64),
                         self.fx, self.fy, self.cx, self.cy, self.alpha)