GeoTIFF / georaster-layer-for-leaflet

Display GeoTIFFs and soon other types of raster on your Leaflet Map
https://geotiff.github.io/georaster-layer-for-leaflet-example/
Apache License 2.0
283 stars 57 forks source link

Geotiff Layer for Anti Meridian South Pacific Leaflet map. #112

Open maritimeWeatherUser opened 1 year ago

maritimeWeatherUser commented 1 year ago

Hi guys i have a test geotiff file that i have loaded as a layer into Leaflet map. The Geotiff has been georeferenced within QGIS using the EPSG:4326 CRS (gdaltranslate and warp) functions and raw geotiff extracted from a GFS forecast netcdf file band (Thredds = Unidata).

The geotiff produced loads onto the Leaflet map BUT on the West of the anti-meridian geotiff layer shows as expected, and on the east of the anti-meridian the layer is missing this half (not rendering). The geotiff layer only shows on the west of the anti-meridian line.

I have set the map CRS to 4326 (or so i think).

Here is the code i use for jscript:

var map = L.map('map').setView([0, 0], 5);

  // add OpenStreetMap basemap
  L.tileLayer('http://{s}.tile.osm.org/{z}/{x}/{y}.png', {
      crs: L.CRS.EPSG4326, 
      attribution: '&copy; <a href="http://osm.org/copyright">OpenStreetMap</a> contributors'
  }).addTo(map);

  // originally from https://globalwindatlas.info/downloads/gis-files
  var url_to_geotiff_file = "https://www.cassavaweather.com/leaflet-kmz-master/NightTest4326.tif";

  fetch(url_to_geotiff_file)
    .then(response => response.arrayBuffer())
    .then(arrayBuffer => {
      parseGeoraster(arrayBuffer).then(georaster => {
        const min = georaster.mins[0];
        const max = georaster.maxs[0];
        const range = georaster.ranges[0];

        // available color scales can be found by running console.log(chroma.brewer);
        console.log(chroma.brewer);
        var scale = chroma.scale("Viridis");

        var layer = new GeoRasterLayer({
            georaster: georaster,
            opacity: 0.7,
            pixelValuesToColorFn: function(pixelValues) {
              var pixelValue = pixelValues[0]; // there's just one band in this raster

              // if there's zero wind, don't return a color
              if (pixelValue === 0) return null;

              // scale to 0 - 1 used by chroma
              var scaledPixelValue = (pixelValue - min) / range;

              var color = scale(scaledPixelValue).hex();

              return color;
            },
            resolution: 256
        });
        console.log("layer:", layer);
        layer.addTo(map);

        map.fitBounds(layer.getBounds());
      });
    });

Here is the output of gdalinfo for the rendered geotif file:

/leaflet-kmz-master$ gdalinfo NightTest4326.tif

Driver: GTiff/GeoTIFF Files: NightTest4326.tif Size is 145, 50 Coordinate System is: GEOGCRS["WGS 84", ENSEMBLE["World Geodetic System 1984 ensemble", MEMBER["World Geodetic System 1984 (Transit)"], MEMBER["World Geodetic System 1984 (G730)"], MEMBER["World Geodetic System 1984 (G873)"], MEMBER["World Geodetic System 1984 (G1150)"], MEMBER["World Geodetic System 1984 (G1674)"], MEMBER["World Geodetic System 1984 (G1762)"], MEMBER["World Geodetic System 1984 (G2139)"], ELLIPSOID["WGS 84",6378137,298.257223563, LENGTHUNIT["metre",1]], ENSEMBLEACCURACY[2.0]], PRIMEM["Greenwich",0, ANGLEUNIT["degree",0.0174532925199433]], CS[ellipsoidal,2], AXIS["geodetic latitude (Lat)",north, ORDER[1], ANGLEUNIT["degree",0.0174532925199433]], AXIS["geodetic longitude (Lon)",east, ORDER[2], ANGLEUNIT["degree",0.0174532925199433]], USAGE[ SCOPE["Horizontal component of 3D system."], AREA["World."], BBOX[-90,-180,90,180]], ID["EPSG",4326]] Data axis to CRS axis mapping: 2,1 Origin = (110.000000000000000,-5.000000000000000) Pixel Size = (1.000000000000000,-1.000000000000000) Metadata: AREA_OR_POINT=Area Image Structure Metadata: INTERLEAVE=PIXEL Corner Coordinates: Upper Left ( 110.0000000, -5.0000000) (110d 0' 0.00"E, 5d 0' 0.00"S) Lower Left ( 110.0000000, -55.0000000) (110d 0' 0.00"E, 55d 0' 0.00"S) Upper Right ( 255.000, -5.000) (255d 0' 0.00"E, 5d 0' 0.00"S) Lower Right ( 255.000, -55.000) (255d 0' 0.00"E, 55d 0' 0.00"S) Center ( 182.500, -30.000) (182d30' 0.00"E, 30d 0' 0.00"S) Band 1 Block=145x14 Type=Byte, ColorInterp=Red Mask Flags: PER_DATASET ALPHA Band 2 Block=145x14 Type=Byte, ColorInterp=Green Mask Flags: PER_DATASET ALPHA Band 3 Block=145x14 Type=Byte, ColorInterp=Blue Mask Flags: PER_DATASET ALPHA Band 4 Block=145x14 Type=Byte, ColorInterp=Alpha

Here is a photo of the results of much trial and error the geotiff layer should be centered on the south pacific region with the anti-meridian at the center of the geotiff layer (achieved). Only the west of the layer from antimeridian shows.

Have tried Stuarts Leaflet plugin but complains about missing dependencies and i can see the files are dated 5 years ago.

Screenshot from 2022-11-11 14-34-31

DanielJDufour commented 1 year ago

Hello, @maritimeWeatherUser . Would you be able to share the GeoTIFF file? It will make debugging easier. If you can't share it publicly on Github, you could email it to daniel.j.dufour@gmail.com. Thank you for bringing this bug to our attention!

maritimeWeatherUser commented 1 year ago

Hi Daniel i have solved the problem through geoserver and openlayers am excited about the future of the project now. As the system i am building is based in the South Pacific one needed to use a specialised CRS system. The QGIS geotiffs where produced correctly. Thanks guys.

On Fri, 11 Nov 2022 at 14:50, Daniel J. Dufour @.***> wrote:

Hello, @maritimeWeatherUser https://github.com/maritimeWeatherUser . Would you be able to share the GeoTIFF file? It will make debugging easier. If you can't share it publicly on Github, you could email it to @.*** Thank you for bringing this bug to our attention!

— Reply to this email directly, view it on GitHub https://github.com/GeoTIFF/georaster-layer-for-leaflet/issues/112#issuecomment-1311127802, or unsubscribe https://github.com/notifications/unsubscribe-auth/AXFTVQKV6U6LZZFOWR6JP5DWHWQ6RANCNFSM6AAAAAAR5EXVAA . You are receiving this because you were mentioned.Message ID: @.***>