Spilhaus projection

[AsgerPetersen]

Also see A story about the "Spilhaus projection" – a map projection that went viral in fall 2018 and will be supported in the next release of ArcGIS

The article finds that the spilhaus projection is a Adams World in a Square II projection with these params:

    Geographic Coordinate System: GCS WGS 1984
    Projection: Adams Square II
    False Easting: 0 m
    False Northing: 0 m
    Scale Factor: 1
    Azimuth: 40.17823482°
    Longitude of Center: 66.94970198°E
    Latitude of Center: 49.56371678°S
    XY Plane Rotation: 45°

It would be fun to have this projection in proj. This document suggests that the "Adams World in a Square II" projection was implemented in the libproject distribution of proj. So maybe it is even possible to find working code somewhere.

[kbevers]

Code is here: https://github.com/jeffbaumes/jeffbaumes-vtk/blob/master/Utilities/vtklibproj4/proj_guyou.c

Have fun translating it to modern day PROJ - Does seem to be five for one kind of deal, so there's that :-)

[AsgerPetersen]

I researched it a bit and wanted to document my findings in case someone would need it later. It seems to be a hot projection at the moment. At least ESRI is making a lot of noise about it at the moment because the new ArcGIS Pro supports it. See for instance this tweet by John Nelson. It even got a dedicated article in forbes a week ago.

[rouault]

Porting Adams Square II from libproj to PROJ wouldn't be hard. But there's more work to do since according to the ESRI story, Spilhaus adds "The Azimuth parameter is a direction from the North towards the top "pole" at the center on the conformal sphere" There would be also other additional work because libproj Adams Square II is for spherical only: "Adams only presented the forward equations for spherical Earth models. Most of today’s geospatial data is defined based on ellipsoidal models, such as WGS 1984 or GRS 1980. That data needs to be converted back from the projected surface to geographic coordinates and we needed to derive the forward and inverse equations for ellipsoidal Earth models. We achieved that by converting geodetic coordinates to a conformal sphere, conformally shrinking the model to a hemisphere, and resolving a complex elliptic integral of the first kind. " So some cartographic expertise required here...

@beuan I see you mentionned as on of the author of the ESRI story. Any plan to contribute it to PROJ ?

[kbevers]

So some cartographic expertise required here...

I am fairly certain that expertise is reachable on our mailing list. You may not get lines of code but I am sure that someone can provide some insigt into as to how to solve the problem.

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[rouault]

Note: if support for Spilhaus is implemented, an action item would be to make sure that ESRI:54099 "WGS_1984_Spilhaus_Ocean_Map_in_Square" which was "blacklisted" in the mapping of ESRI WKT Adams_Square_II in https://github.com/OSGeo/PROJ/pull/2157 is enabled

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[philippemiron]

I'm not familiar enough with PROJ, but I found this implementation. Is this something that could work if translated to proj syntax?

Adams World in a Square I + II

Porting Adams Square II from libproj to PROJ wouldn't be hard. But there's more work to do since according to the ESRI story, Spilhaus adds "The Azimuth parameter is a direction from the North towards the top "pole" at the center on the conformal sphere" There would be also other additional work because libproj Adams Square II is for spherical only: "Adams only presented the forward equations for spherical Earth models. Most of today’s geospatial data is defined based on ellipsoidal models, such as WGS 1984 or GRS 1980. That data needs to be converted back from the projected surface to geographic coordinates and we needed to derive the forward and inverse equations for ellipsoidal Earth models. We achieved that by converting geodetic coordinates to a conformal sphere, conformally shrinking the model to a hemisphere, and resolving a complex elliptic integral of the first kind. " So some cartographic expertise required here...

@beuan I see you mentionned as on of the author of the ESRI story. Any plan to contribute it to PROJ ?

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[AJueling]

@rouault, this question may reveal my ignorance, but would it not be good enough to have the projection work only from spherical data for the majority of use cases which would be oceanographic and global in scope? I think this projection would be of great interest to oceanographers for predominantly illustrative (and perhaps less scientific) purposes.

[rouault]

would it not be good enough to have the projection work only from spherical data for the majority of use cases which would be oceanographic and global in scope?

A spherical version would be a good start, but someone would have to find the maths needed to implement the Azimuth and Latitude of Center parameters needed for Spilhaus, and not present in the PROJ implementation of Adams World in a Square II

[kerkovits]

would it not be good enough to have the projection work only from spherical data for the majority of use cases which would be oceanographic and global in scope?

A spherical version would be a good start, but someone would have to find the maths needed to implement the Azimuth and Latitude of Center parameters needed for Spilhaus, and not present in the PROJ implementation of Adams World in a Square II

This should be extremely easy. The map projection is conformal. If you want to map ellipsoidal data just use conformal latitude (formula 3-1 of https://pubs.usgs.gov/pp/1395/report.pdf) instead of geographic latitude and just use spherical formulae. This will be 100% conformal and works with every conformal map projection.

The "center of projection" and "azimuth" is just a graticule rotation before applying the map projection. (formula 5-7 and 5-8b of https://pubs.usgs.gov/pp/1395/report.pdf) I could reproduce the projection with alpha=30°, lambda_0=115°. Unfortunately the exact value of beta is unknown (it should be calculated using spherical trigonometry) but I could achieve an almost perfect result with beta=-28.8° (determined by trial and error).

So this should be possible.

[philippemiron]

would it not be good enough to have the projection work only from spherical data for the majority of use cases which would be oceanographic and global in scope?

A spherical version would be a good start, but someone would have to find the maths needed to implement the Azimuth and Latitude of Center parameters needed for Spilhaus, and not present in the PROJ implementation of Adams World in a Square II

This should be extremely easy. The map projection is conformal. If you want to map ellipsoidal data just use conformal latitude (formula 3-1 of https://pubs.usgs.gov/pp/1395/report.pdf) instead of geographic latitude and just use spherical formulae. This will be 100% conformal and works with every conformal map projection.

The "center of projection" and "azimuth" is just a graticule rotation before applying the map projection. (formula 5-7 and 5-8b of https://pubs.usgs.gov/pp/1395/report.pdf) I could reproduce the projection with alpha=30°, lambda_0=115°. Unfortunately the exact value of beta is unknown (it should be calculated using spherical trigonometry) but I could achieve an almost perfect result with beta=-28.8° (determined by trial and error).

So this should be possible.

Do you have a a working example?

[kerkovits]

Do you have a a working example?

Yes, I do but it is outside the PROJ ecosystem, and is available only in Hungarian. This program was created for personal use, so do not expect easy-to-understand program codes. I can only share the mathematical formula, I am bad at programming.

http://mercator.elte.hu/~kerkovits/projections/megnez.php?vetulet=adams&s1=115&s2=30&s3=-28.8

The image is rotated by 45° but it is not a big issue.

[kerkovits]

I found that the parameter beta in the formula I linked before should be atan2(-sin(azimuth),-tan(latitude_of_center)) This gives beta=-28.8013495219182°, which is fairly close to the value determined by trial and error.

[philippemiron]

Can we reopened this issue? And is there someone that could help me with this?

[kbevers]

@philippemiron if you are willing to take this on I can help you along as time allow it

[kerkovits]

I did some calculations and it seems that ellipsoidal formulae are a bit more complicated than I originally assumed but still possible to implement:

Step 1: calculate conformal latitude of of the center. Use formula (3-1) of https://pubs.usgs.gov/pp/1395/report.pdf substitute the latitude of center for phi, and the first eccentricity to e. I will call the result chi as the conformal_latitude_of_center. This step should be omitted for spherical calculations.

Step 2: calculate the following quantities: alpha=asin(cos(conformal_latitude_of_center)*cos(azimuth)) (should be a bit less than +30° for ellipsoidal formulae and exactly +30° for a sphere) lambda_0=longitude_of_center+atan2(tan(azimuth),-sin(conformal_latitude_of_center)) (should be exactly +115°) beta=atan2(-sin(azimuth),-tan(conformal_latitude_of_center)) (should be approximately -29°)

Step 3: Calculate the conformal latitude of the mapped point by formula (3-1). Skip this step for spherical formula.

Step 4: Calculate phi_prime and lambda_prime from formulae (5-7) and (5-8b). Use constants alpha, beta, and lambda_0 from step 2 and use the conformal latitude from step 3 for phi. Lambda is the longitude of the mapped point.

Step 5: Apply the spherical Adams projection for phi_prime and lambda_prime.

Step 6: Rotate the result by 45° clockwise.

Step 7: Rejoice.

[kerkovits]

P.S. You should use atan2 for the arctan in formula (5-8b) at step 4.

[iapafoto]

Java implementation (Proj4J) Thanks to kerkovits comments on https://github.com/OSGeo/PROJ/issues/1851 Formulas ref [Map projections - a working Manual] by John P.Snyder

private final static double
        α = Math.toRadians(30), sinα = sin(α), cosα = cos(α),
        β = Math.toRadians(-28.8013495219182);

public ProjCoordinate changePole(double λ, double ϕ, ProjCoordinate out) {
    final double 
            cosϕ = cos(-ϕ), sinϕ = sin(-ϕ),
            cosλ = cos(-λ), sinλ = sin(-λ),
            sin_ϕp = sinα*sinϕ - cosα*cosϕ*cosλ, //  5-7
            ϕp = asin(sin_ϕp),
            λp = atan2(cosϕ*sinλ, sinα*cosϕ*cosλ + cosα*sinϕ) + β; // 5-8b
    out.x = -ProjectionMath.normalizeLongitude(λp);
    out.y = -ϕp;
    return out;
}

@Override
public ProjCoordinate project(double lplam, double lpphi, ProjCoordinate xy) {
    // Conformal latitude
    double x = 2.*atan(tan(Math.PI/4. + lpphi/2)*pow((1 - e*sin(lpphi))/(1 + e*sin(lpphi)),e/2))- Math.PI/2.; // 3-1
    final ProjCoordinate rotLamPhi = changePole(lplam, x, new ProjCoordinate());
    lplam = rotLamPhi.x;
    lpphi = rotLamPhi.y;
    final double 
            spp = tan(.5 * lpphi),
            sa = cos(aasin(spp)) * sin(.5 * lplam),
            b = aacos(spp),
            a = aacos(cos(aasin(spp)) * sin(.5 * lplam));
    xy.x = ell_int_5(((spp + sa) < 0. ? -1 : 1) * aasin(sqrt(   (1. + min(0., cos(a + b))))));
    xy.y = ell_int_5(((spp - sa) < 0. ? -1 : 1) * aasin(sqrt(abs(1. - max(0., cos(a - b))))));
    return xy;
}

[rouault]

Java implementation (Proj4J)

Is it officially integrated ? I can't spot it in https://github.com/locationtech/proj4j/tree/master/src/main/java/org/locationtech/proj4j/proj

[iapafoto]

no I just start implementing it on a local version (a lot of projections are missing or implemented with error in Proj4J)

[rouault]

Could someone monitoring this issue create some test points (longitude, latitude) -> (X,Y) with Esri software, possibly using ESRI:54099 "WGS_1984_Spilhaus_Ocean_Map_in_Square", so we can check interoperability ?

[iapafoto]

For any useful purpose, a full code with projection and inverse projection (Java version) Notice: Reverse projection loses precision when approaching [latitude -25..-30 longitude 175..180] and [latitude 30..35 longitude 0..5]

`/* Copyright 2021 Sebastien Durand based on (2006 Jerry Huxtable)

Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except in compliance with the License. You may obtain a copy of the License at

http://www.apache.org/licenses/LICENSE-2.0

Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License. */

package org.osgeo.proj4j.proj;

import static java.lang.Math.abs; import static java.lang.Math.asin; import static java.lang.Math.atan2; import static java.lang.Math.cos; import static java.lang.Math.max; import static java.lang.Math.min; import static java.lang.Math.sin; import static java.lang.Math.sqrt; import static java.lang.Math.tan; import org.osgeo.proj4j.ProjCoordinate; import org.osgeo.proj4j.parser.Proj4Keyword; import org.osgeo.proj4j.util.ProjectionMath; import static org.osgeo.proj4j.util.ProjectionMath.aacos; import static org.osgeo.proj4j.util.ProjectionMath.aasin;

public class SpilhausProjection extends Projection {

protected final static double TOL = 1e-9;
protected final static double RSQRT2 = 0.7071067811865475244008443620;
protected final static double M_2_PI = 0.63661977236758134308;
protected final static double SEC_TO_RAD = Math.PI/(3600.*180.); 

private final static double C[] = {
    -8.58691003636495e-07,
    2.02692115653689e-07,
    3.12960480765314e-05,
    5.30394739921063e-05,
    -0.0012804644680613,
    -0.00575574836830288,
    0.0914203033408211,};

private static double ell_int_5(double phi) {
    /* Procedure to compute elliptic integral of the first kind
     * where k^2=0.5.  Precision good to better than 1e-7
     * The approximation is performed with an even Chebyshev
     * series, thus the coefficients below are the even values
     * and where series evaluation  must be multiplied by the argument. */
    double C0 = 2.19174570831038;
    double y = phi * M_2_PI;
    y = 2. * y * y - 1.;
    double y2 = 2. * y, d1 = 0., d2 = 0.0;
    for (double c : C) {
        double temp = d1;
        d1 = y2 * d1 - d2 + c;
        d2 = temp;
    }
    return phi * (y * d1 - d2 + 0.5 * C0);
}

public SpilhausProjection() {
    initialize();
}

/**
 * Applying these relationships to transformations, without showing some intermediate derivations, 
 * formulas (5-7) through (5-8b) are obtained. 
 * To place the North Pole of the sphere at a latitude α on a meridian β east of the central meridian (λp = 0) 
 * of the basic projection (see fig. 6), the transformed latitude ϕp and transformed longitude λp 
 * on the basic projection which correspond to latitude ϕ and longitude λ of the spherical Earth may be calculated as follows, 
 * letting the central meridian λ0 correspond with λp = β
*/

// Formulas ref from (Map projections - a working Manual) by John P.Snyder
// Thanks to kerkovits comments on https://github.com/OSGeo/PROJ/issues/1851
protected final static double
        α = Math.toRadians(30), sinα = sin(α), cosα = cos(α),
        β = Math.toRadians(-28.8013495219182);

public ProjCoordinate changePole(double λ, double ϕ, ProjCoordinate out) {
    final double 
            λ0 = 0,
            cosϕ = cos(ϕ), sinϕ = sin(ϕ),
            cosλ_λ0 = cos(λ-λ0), sinλ_λ0 = sin(λ-λ0),
            sinϕp = sinα*sinϕ - cosα*cosϕ*cosλ_λ0, //  5-7
            ϕp = asin(sinϕp),
            λp = atan2(cosϕ*sinλ_λ0, sinα*cosϕ*cosλ_λ0 + cosα*sinϕ) + β; // 5-8b
    out.x = ProjectionMath.normalizeLongitude(λp);
    out.y = ϕp;
    return out;
}

public ProjCoordinate changePoleInv(double λp, double ϕp, ProjCoordinate out) {
    final double
            λ0 = 0,
            cosλp_β = cos(λp-β),
            cosϕp = cos(ϕp), sinϕp = sin(ϕp),
            ϕ = asin(sinα*sinϕp + cosα*cosϕp*cosλp_β),                        // 5-9
            λ = atan2(cosϕp*sin(λp-β),(sinα*cosϕp*cosλp_β-cosα*sinϕp)) + λ0;  // 5-10b
    out.x = ProjectionMath.normalizeLongitude(λ);  
    out.y = ϕ;
    return out;
}

public static ProjCoordinate changePoleFull(double λ, double ϕ, double α, double β, double λ0, ProjCoordinate out) {
    final double 
            cosϕ = cos(ϕ), sinϕ = sin(ϕ),
            cosλ = cos(λ-λ0), sinλ = sin(λ-λ0),
            sin_ϕp = sin(α)*sinϕ - cos(α)*cosϕ*cosλ, //  5-7
            ϕp = asin(sin_ϕp),
            λp = atan2(cosϕ*sinλ, sin(α)*cosϕ*cosλ + cos(α)*sinϕ) + β; // 5-8b
    out.x = ProjectionMath.normalizeLongitude(λp);
    out.y = -ϕp;
    return out;
}

@Override
public ProjCoordinate project(double lplam, double lpphi, ProjCoordinate xy) {
    // Conformal latitude
    double x = lpphi - SEC_TO_RAD*(700.0427 * sin(2 * lpphi) + 0.9900 * sin(4 * lpphi) + 0.0017 * sin(6 * lpphi));  // 3-3
    // Change pole
    final ProjCoordinate rotLamPhi = changePole(-lplam, -x, new ProjCoordinate());
    // Admas Square 2 Projection
    return projectAdamsSq2(-rotLamPhi.x, -rotLamPhi.y, xy);
}        

public static ProjCoordinate projectAdamsSq2(double lplam, double lpphi, ProjCoordinate xy) {
    final double 
            spp = tan(.5 * lpphi),
            sa = cos(aasin(spp)) * sin(.5 * lplam),
            b = aacos(spp),
            a = aacos(cos(aasin(spp)) * sin(.5 * lplam));
    xy.x = ell_int_5(((spp + sa) < 0. ? -1 : 1) * aasin(sqrt(   (1. + min(0., cos(a + b))))));
    xy.y = ell_int_5(((spp - sa) < 0. ? -1 : 1) * aasin(sqrt(abs(1. - max(0., cos(a - b))))));
    return xy;
}

@Override
public ProjCoordinate projectInverse(double xyx, double xyy, ProjCoordinate out) {
    // Inverse Adams square II projection
    out.y = max(min(xyy / 2.62181347, 1.0), -1.0) * M_HALFPI;
    out.x = abs(out.y) >= M_HALFPI ? 0 : max(min(xyx / 2.62205760 / cos(out.y), 1.0), -1.0) * M_PI;
    ProjCoordinate out2 = pj_generic_inverse_2d_AdamsSq2(xyx, xyy, out);
    // Inverse change Pole
    changePoleInv(-out2.x, -out2.y, out);
    // Inverse Conformal 
    double x = -out.y;
    out.y = x + SEC_TO_RAD*(700.0420* sin(2*x) + 1.3859 *sin(4*x) + 0.0037 * sin(6*x));   // 3-6
    out.x = -out.x;
    return out;
}

@Override
public void initialize() {
    super.initialize();
    this.es = 0;
    this.projectionLongitude = Math.toRadians(115);
}

/**
 * Returns true if this projection is equal area
 */
@Override
public boolean isEqualArea() {
    return true;
}

@Override
public boolean hasInverse() {
    return true;
}

/**
 * Returns the ESPG code for this projection, or 0 if unknown.
 *
 * @return
 */
@Override
public int getEPSGCode() {
    return 54099;
}

@Override
public String toString() {
    return "Spilhaus";
}

@Override
public Proj4Keyword[] getLstRequiredParams() {
    return new Proj4Keyword[]{};
}

@Override
public Proj4Keyword[] getLstOptinalParams() {
    return new Proj4Keyword[]{};
}

public ProjCoordinate pj_generic_inverse_2d_AdamsSq2(double xyx, double xyy, ProjCoordinate lpInitial) {
    ProjCoordinate lp = new ProjCoordinate(lpInitial.x, lpInitial.y);
    double deriv_lam_X = 0;
    double deriv_lam_Y = 0;
    double deriv_phi_X = 0;
    double deriv_phi_Y = 0;

    for (int i = 0; i < 15; i++) {
        ProjCoordinate xyApprox = projectAdamsSq2(lp.x, lp.y, new ProjCoordinate());

        double deltaX = xyApprox.x - xyx;
        double deltaY = xyApprox.y - xyy;
        if (abs(deltaX) < 1e-10 && abs(deltaY) < 1e-10) {
            return lp;
        }

        if (i == 0 || abs(deltaX) > 1e-6 || abs(deltaY) > 1e-6) {
            // Compute Jacobian matrix (only if we aren't close to the final
            // result to speed things a bit)
            ProjCoordinate lp2 = new ProjCoordinate();
            ProjCoordinate xy2 = new ProjCoordinate();
            double dLam = lp.x > 0 ? -1e-6 : 1e-6;
            lp2.x = lp.x + dLam;
            lp2.y = lp.y;
            xy2 = projectAdamsSq2(lp2.x,lp2.y, xy2);
            double deriv_X_lam = (xy2.x - xyApprox.x) / dLam;
            double deriv_Y_lam = (xy2.y - xyApprox.y) / dLam;

            double dPhi = lp.y > 0 ? -1e-6 : 1e-6;
            lp2.x = lp.x;
            lp2.y = lp.y + dPhi;
            xy2 = projectAdamsSq2(lp2.x,lp2.y, xy2);
            double deriv_X_phi = (xy2.x - xyApprox.x) / dPhi;
            double deriv_Y_phi = (xy2.y - xyApprox.y) / dPhi;

            // Inverse of Jacobian matrix
            double det = deriv_X_lam * deriv_Y_phi - deriv_X_phi * deriv_Y_lam;
            if (det != 0) {
                deriv_lam_X = deriv_Y_phi / det;
                deriv_lam_Y = -deriv_X_phi / det;
                deriv_phi_X = -deriv_Y_lam / det;
                deriv_phi_Y = deriv_X_lam / det;
            }
        }

        if (xyx != 0) {
            // Limit the amplitude of correction to avoid overshoots due to
            // bad initial guess
            double delta_lam = max( min(deltaX * deriv_lam_X + deltaY * deriv_lam_Y, 0.3), -0.3);
            lp.x -= delta_lam;
            if (lp.x < -M_PI)
                lp.x = -M_PI;
            else if (lp.x > M_PI)
                lp.x = M_PI;
        }

        if (xyy != 0) {
            double delta_phi = max( min(deltaX * deriv_phi_X + deltaY * deriv_phi_Y, 0.3), -0.3);
            lp.y -= delta_phi;
            if (lp.y < -M_HALFPI)
                lp.y = -M_HALFPI;
            else if (lp.y > M_HALFPI)
                lp.y = M_HALFPI;
        }
    }
  //  proj_context_errno_set(P->ctx, PROJ_ERR_COORD_TRANSFM_OUTSIDE_PROJECTION_DOMAIN);
    return lp;
}

} `

[rouault]

For any useful purpose, a full code with projection and inverse projection (Java version)

@iapafoto : is it your code ? (I'm confused by the "Copyright 2006 Jerry Huxtable") Could that be used under the X/MIT license used by PROJ ? Some part (the inversion part) seems to come from PROJ.

[iapafoto]

Yes it is my code as an addon to Proj4J, but it is based on Proj4 C version of Adams Square II projection (https://github.com/OSGeo/PROJ/blob/master/src/projections/adams.cpp)

The Copyright was a copy-paste of other Proj4J projections

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[nyalldawson]

For reference, here's a js implementation: https://github.com/neocarto/bertin/blob/main/src/projections/spilhaus.js

[mdsumner]

Could someone monitoring this issue create some test points (longitude, latitude) -> (X,Y) with Esri software, possibly using ESRI:54099 "WGS_1984_Spilhaus_Ocean_Map_in_Square", so we can check interoperability ?

inadvertently, not seeing this comment I did this here: https://github.com/mdsumner/spilheist/tree/main/inst/extdata

(Matthew Law provided the files, I wanted fields of lon,lat,x,y to do empirical lookup)

[fabeit]

Can someone please tell me how to use this projection in R? I have tried to take the WKT and WKT2 and substitute into a crs object but I get this error "Scale limits cannot be mapped onto spatial coordinates in coord_sf()."

[kbevers]

Can someone please tell me how to use this projection in R?

Get in touch with whoever has implemented it in R

[fabeit]

@kbevers I can't find any information on the implementation in R, there is an old question on GISexchange which is not answered. I can read the projection details correctly but the plotting doesn't work.

[kbevers]

can't find any information on the implementation in R

You are not very likely to find it here either

[snowman2]

I recommend going here with R questions related to PROJ:

• https://github.com/hypertidy/PROJ
• https://cran.r-project.org/web/packages/PROJ/index.html
• https://cran.r-project.org/web/packages/proj4/proj4.pdf

[rtlemos]

Here's an implementation in R and Python: https://github.com/rtlemos/spilhaus/tree/main

[noirchen]

@rtlemos This helps a lot.