Contour Layer crashes on coastline

[gsueur]

Hi,

I've set up a contour layer using a VRT file assembling 4 Copernicus EU DEM files.
The contour layer is neat, but it makes MapServer to crash every time I go on the coastline.
How can I get a complete error message to help identify the problem ?

Thanks

[rouault]

Run the request with
MS_ERRORFILE=log.txt MS_DEBUGLEVEL=5 [valgrind | gdb --args] mapserv QUERY_STRING="..."

[gsueur]

Actually it just freezes like this:
[Sat Nov 21 15:34:29 2020].204032 msLoadMap(): 0.015s
[Sat Nov 21 15:34:29 2020].206380 msDrawMap(): rendering using outputformat named png (AGG/PNG).
[Sat Nov 21 15:34:29 2020].206417 Entering msContourLayerOpen().
[Sat Nov 21 15:34:29 2020].208653 Entering msContourLayerReadRaster().
[Sat Nov 21 15:34:29 2020].208679 msContourLayerReadRaster(): Entering transform.
[Sat Nov 21 15:34:29 2020].209165 msContourLayerReadRaster(): src=24977,39552,460,332, dst=0,0,460,332

[jmckenna]

@gsueur I can reproduce here ("Killed") with a single Copernicus file with GDAL 2.4.4, and with GDAL 3.2.0 it at least gives a better message ("GDALContourGenerate() failed: Input values and/or interval settings would lead to too many levels.") I also set inside my mapfile CONFIG "CPL_DEBUG" "ON"

This error happens if you only set a single CONTOUR_INTERVAL, and not handling outer scales.

Try again but replace your PROCESSING "CONTOUR_INTERVAL= ... with these 3 lines:

    PROCESSING "CONTOUR_INTERVAL=0,25000:5" # interval of 5 for scales of 25000 or less
    PROCESSING "CONTOUR_INTERVAL=25000,500000:10" # interval of 10 for scales in the 25000 to 500000 range
    PROCESSING "CONTOUR_LEVELS=500000,0:10,25,50,100,200,300,400,500,600,700,800,900,1000" # explicit list of levels for scales of 500000 and up 

[jmckenna]

For the record, here is my working test LAYER (works with both GDAL 2.4.4 and GDAL 3.2.0) :

  /* Contours (generated from Copernicus EU DEM file) */
  LAYER
    NAME "dem-contours"
    TYPE LINE
    STATUS ON
    CONNECTIONTYPE CONTOUR
    DATA "eu_dem_v11_E40N20_clip.TIF"
    PROCESSING "BANDS=1"
    PROCESSING "CONTOUR_ITEM=Band_1"
    PROCESSING "CONTOUR_INTERVAL=0,25000:5" # interval of 5 for scales of 25000 or less
    PROCESSING "CONTOUR_INTERVAL=25000,500000:10" # interval of 10 for scales in the 25000 to 500000 range
    PROCESSING "CONTOUR_LEVELS=500000,0:10,25,50,100,200,300,400,500,600,700,800,900,1000" # explicit list of levels for scales of 500000 and up 
    CLASS
      NAME "Contours"
      STYLE
        WIDTH 1
        COLOR 139 69 19
      END
    END
  END # layer

Attached coastline (Corsica island) contours from MapServer, and then zoomed in to inside island...

[gsueur]

Hi Jeff,
Thanks for giving a head. Actually your trick didn't make it here. I feel a bit stupid not having posted my LAYER before, here it is

 LAYER
    GROUP "default"
    NAME "contours14"
    TYPE LINE
    STATUS ON
    PROJECTION
        "init=epsg:3035"
    END
    MAXSCALEDENOM 40626
    MINSCALEDENOM 20313
     METADATA
        "wms_extent" "2426378.0132, 1528101.2618, 6293974.6215, 5446513.5222"
     END
    CONNECTIONTYPE CONTOUR
    DATA "/opt/neogeo/dem/aster_dem.vrt"
    PROCESSING "BANDS=1"
    PROCESSING "CONTOUR_ITEM=elev"
    PROCESSING "CONTOUR_INTERVAL=0,25000:10"
    PROCESSING "CONTOUR_INTERVAL=25000,100000:20"
    PROCESSING "CONTOUR_LEVELS=100000,0:25,50,100,200,300,400,500,600,700,800,900,1000"
    GEOMTRANSFORM (smoothsia(generalize([shape], 0.25*[data_cellsize])))
    LABELITEM "elev"
    CLASS
        EXPRESSION ([elev] % 50 = 0)
        STYLE
            COLOR 136 99 85
            WIDTH 0.8
        END
        LABEL
          TYPE TRUETYPE
          PARTIALS FALSE
          FONT nsli
          COLOR 136 99 85
          OUTLINECOLOR 255 255 255
          OUTLINEWIDTH 2
          TEXT (tostring([elev],"%.0f"))
          SIZE 7
          ANGLE FOLLOW MAXOVERLAPANGLE 25
          TYPE TRUETYPE
          MAXOVERLAPANGLE 22.5
          REPEATDISTANCE 400
          MINDISTANCE 100
          OFFSET 0 10
        END
    END
    CLASS
        EXPRESSION ([elev] > 0)
        STYLE
            COLOR 136 99 85
            WIDTH 0.4
        END
    END
END

but none of the instructions (even the geomtransform thing) seems to be the root of the problem.
Could it come from NODATA values not properly cleaned in my DEM ?
Best Regards

[jmckenna]

@gsueur here is a small snippet of your Copernicus EU data, with mapfile, can you try this locally on your machine? (see 'commands.txt' file for the shp2img call. Also see line#10 of the mapfile, to enable that coastline EXTENT). After you get a valid image, then add your missing additions from your other mapfile, I bet you can find what line crashes that way. https://gatewaygeomatics.com/dl/ticket-6182-contour.zip

[jmckenna]

@gsueur but indeed, next time you can save me all that time to download the 5GB dem file, clip it with gdalwarp, create new mapfile, create commands.txt file with commands to reproduce problem.........by doing that all yourself and attaching the small .zip here to the ticket. Cheers :)

[gsueur]

Sorry Jeff, I didn't mean to bother you. All I was asking for was some command to help me investigate, not you to completely build a test case for me.
I have uploaded and tested it. I have the pretty same result : correct display on inlands and/or small scales. No image and a hanging process on coastlines. Apparently, with
CONTOUR_INTERVAL=25000,500000:10
CONTOUR_LEVELS=500000,0:...
it breaks precisely at scale 1/1000000
but with
CONTOUR_INTERVAL=25000,50000:10
CONTOUR_LEVELS=50000,0:...
it breaks at scale 1/100000 !!
So we have got something here !
Digging into more complete investigations, will report my findings if any.
Have a great day and please don't be upset.
Best regards
Guillaume

[gsueur]

what I have so far :

• there is a factor 2 effect between the requested scale in CONTOUR_LEVELS and the actual scale which makes it fire in QGIS (using EPSG:3035 to avoid projection noise). It can be seen with CONTOUR_LEVELS=5000,0:10,20,30,40 --> you have 4 contours lines util you cross the 1/10000 scale into higher scale, where many more lines are then drawn (NO CONTOUR_INTERVAL used in the mapfile). This factor effect is not a big deal, but you need to have it in mind to understand the next point.
• Once you get into any CONTOUR_INTERVAL range on coastlines, it hangs.
• If I use a single and simple PROCESSING "CONTOUR_INTERVAL=10", it hangs.

Can't we have a division by zero here somewhere ?

[gsueur]

Valgrind output :
==16233== Using Valgrind-3.14.0 and LibVEX; rerun with -h for copyright info
==16233== Command: mapserv QUERY_STRING=map=/maps/contour.map&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=2158113.916165928356,4267809.092852747999,2168384.693913876079,4294272.538221334107&CRS=EPSG:3035&WIDTH=1430&HEIGHT=555&LAYERS=dem-contours&STYLES=&FORMAT=image/png&DPI=72&MAP_RESOLUTION=72&FORMAT_OPTIONS=dpi:72&TRANSPARENT=TRUE
==16233==
vex amd64->IR: unhandled instruction bytes: 0x62 0xF1 0xFE 0x8 0x6F 0x2 0xC5 0xF8 0x11 0x44
vex amd64->IR: REX=0 REX.W=0 REX.R=0 REX.X=0 REX.B=0
vex amd64->IR: VEX=0 VEX.L=0 VEX.nVVVV=0x0 ESC=NONE
vex amd64->IR: PFX.66=0 PFX.F2=0 PFX.F3=0
==16233== valgrind: Unrecognised instruction at address 0x48d0b78.
==16233== at 0x48D0B78: msIO_vfprintf (in /usr/lib/libmapserver.so.7.6.1)
==16233== by 0x48D0DBE: msIO_fprintf (in /usr/lib/libmapserver.so.7.6.1)
==16233== by 0x49079DB: msDebug (in /usr/lib/libmapserver.so.7.6.1)
==16233== by 0x48867B5: loadParams (in /usr/lib/libmapserver.so.7.6.1)
==16233== by 0x1092EF: main (in /usr/bin/mapserv)
==16233== Your program just tried to execute an instruction that Valgrind
==16233== did not recognise. There are two possible reasons for this.
==16233== 1. Your program has a bug and erroneously jumped to a non-code
==16233== location. If you are running Memcheck and you just saw a
==16233== warning about a bad jump, it's probably your program's fault.
==16233== 2. The instruction is legitimate but Valgrind doesn't handle it,
==16233== i.e. it's Valgrind's fault. If you think this is the case or
==16233== you are not sure, please let us know and we'll try to fix it.
==16233== Either way, Valgrind will now raise a SIGILL signal which will
==16233== probably kill your program.
==16233==
==16233== Process terminating with default action of signal 4 (SIGILL)
==16233== Illegal opcode at address 0x48D0B78
==16233== at 0x48D0B78: msIO_vfprintf (in /usr/lib/libmapserver.so.7.6.1)
==16233== by 0x48D0DBE: msIO_fprintf (in /usr/lib/libmapserver.so.7.6.1)
==16233== by 0x49079DB: msDebug (in /usr/lib/libmapserver.so.7.6.1)
==16233== by 0x48867B5: loadParams (in /usr/lib/libmapserver.so.7.6.1)
==16233== by 0x1092EF: main (in /usr/bin/mapserv)
==16233==
==16233== HEAP SUMMARY:
==16233== in use at exit: 264,053 bytes in 843 blocks
==16233== total heap usage: 1,457 allocs, 614 frees, 361,969 bytes allocated
==16233==
==16233== LEAK SUMMARY:
==16233== definitely lost: 0 bytes in 0 blocks
==16233== indirectly lost: 0 bytes in 0 blocks
==16233== possibly lost: 0 bytes in 0 blocks
==16233== still reachable: 264,053 bytes in 843 blocks
==16233== suppressed: 0 bytes in 0 blocks
==16233== Rerun with --leak-check=full to see details of leaked memory
==16233==
==16233== For counts of detected and suppressed errors, rerun with: -v
==16233== ERROR SUMMARY: 0 errors from 0 contexts (suppressed: 0 from 0)
Illegal instruction

[gsueur]

gdb output :
Starting program: /usr/lib/cgi-bin/mapserv QUERY_STRING="map=/maps/contour.map&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=2158113.916165928356,4267809.092852747999,2168384.693913876079,4294272.538221334107&CRS=EPSG:3035&WIDTH=1430&HEIGHT=555&LAYERS=dem-contours&STYLES=&FORMAT=image/png&DPI=72&MAP_RESOLUTION=72&FORMAT_OPTIONS=dpi:72&TRANSPARENT=TRUE"
[Thread debugging using libthread_db enabled]
Using host libthread_db library "/lib/x86_64-linux-gnu/libthread_db.so.1".
GNM: GNMRegisterAllInternal
GNM: RegisterGNMFile
GNM: RegisterGNMdatabase
GDAL: GDALOpen(/maps/./data/eu_dem_v11_E40N20_clip.TIF, this=0x55555563c160) succeeds as GTiff.
GDAL: GDAL_CACHEMAX = 1591 MB
GDAL: GDALOpen(MEM:::DATAPOINTER=0x7fffee6b8010,PIXELS=1070,LINES=423,BANDS=1,DATATYPE=Float64, this=0x5555571d4e10) succeeds as MEM.
GNM: GNMRegisterAllInternal
GNM: RegisterGNMFile
GNM: RegisterGNMdatabase
GDAL: GDALDriver::Create(Memory,,0,0,0,Unknown,(nil))

And then it hangs forever.

[jmckenna]

@gsueur so my test case crashes for you? Here I cannot reproduce on several machines and GDAL/MapServer versions, and on different operating systems. What EXTENT values are you using in my mapfile (see line#10) ?

You report "If I use a single and simple PROCESSING "CONTOUR_INTERVAL=10", it hangs" <--- please see my feedback above, as your single/simple line will crash MapServer each time, and i explained above how recent GDAL versions handle the issue properly. I also updated the docs with this warning, so please avoid the single/simple line as it will overwhelm MapServer/GDAL/OGR at zoomed-out scales.

I'm not upset, but please do keep in mind when reporting issues to take the time to create a small test case & attach it to your ticket report, it makes it easier for others to test your problems. :) Enjoy your Sunday.

[jmckenna]

ah I can see your extent values used in your output....will test here with the test case but with your EXTENT values, hopefully I see a crash too, ha

[jmckenna]

Can reproduce here. The EXTENT value (around line#10 of mapfile) to reproduce is:

  EXTENT 4267809.092852747999 2158113.916165928356 4294272.538221334107 2168384.693913876079 #Guillaume's

and GDAL 3.2.0 gives the error:

GDALContourGenerate() failed: Input values and/or interval settings would lead to too many levels

These extents are very zoomed-out, for generating contours, I am not surprised there are problems, personally.

Will try to find the correct settings to avoid this error...

[jmckenna]

this is tricky @gsueur !! I see what you mean now about problems around scale 1/100000

I'm trying to handle it nicely through CONTOUR_INTERVAL but I keep getting "Input values and/or interval settings would lead to too many levels"

[jmckenna]

I tend to think this is a bug now, as you said, as I can't seem to handle nicely the outer scales without GDAL/OGR returning the "too many levels" error.

[gsueur]

Indeed, tricky and nasty ! I don't understand why the bigger scales would return "too many levels" when it works fine at lower ones. There is maybe something happening with NODATA causing this, when most of the extent is covering with null values ?

[jmckenna]

very nasty ha. I guess we need to look into the logic in GDALContourGenerate(), to understand this "too many levels" message. But I'm glad you reported this here (!!), as I think whatever we find will help prevent others from hitting this same wall again.

[gsueur]

this ?

https://github.com/OSGeo/gdal/blob/3a5b462edecd6fceb9961e3b3314b4d59c362b46/gdal/alg/marching_squares/level_generator.h

I can't really get into the logic here.

[jmckenna]

I think it would be good also to see how the utility 'gdal_contour' handles this same data.

[gsueur]

I've had generated a whole contour dataset before using the mapserver contour layer. It went ok, as it goes OK when you don't set any CONTOUR_INTERVAL and just rely on CONTOUR_LEVELS. I think there is a computation made for the contour lines based on CONTOUR_INTERVAL options which fires the exception when some pixels have NODATA.

[gsueur]

in mapcontour.c, line 617 :
eErr = GDALContourGenerate( hBand, interval, 0.0,
levelCount, levels,
FALSE, 0.0, hLayer,
OGR_FD_GetFieldIndex(OGR_L_GetLayerDefn( hLayer),
"ID" ),
(elevItem == NULL) ? -1 :
OGR_FD_GetFieldIndex(OGR_L_GetLayerDefn( hLayer),
elevItem ),
NULL, NULL );

why is the option bUseNoData set to FALSE ?

[gsueur]

I think I have a clue : the NODATA value is -3.40282e+38 in the DEM files !!!
that could explain why GDAL complains about the too many lines to generate right ?

[rouault]

I think I have a clue : the NODATA value is -3.40282e+38 in the DEM files !!!
that could explain why GDAL complains about the too many lines to generate right ?

Indeed (if that nodata value is actually hit). Should be fixed per #6183 but I only made sure it compiled. Could you confirm it fixes your issue ?

[gsueur]

Complete crash even on the getCapabilties... sorry :-(

[rouault]

Complete crash even on the getCapabilties... sorry :-(

I'd suspect an issue in your build environment. msautotest is happy with my change and I don't think the modified code paths should be hit by GetCapabilities

[gsueur]

indeed, my server is broken, I can't figure out why. I will fix that tomorrow. Cheers.

[gsueur]

Mapserver fixed, but contour lines still don't show up on the coast.

[gsueur]

Actually the error is correctly handled by MapServer now, which doesn't crash no more.
But the calculation of the number of levels must still be erroneous.:
:


msDrawMap(): Image handling error. Failed to draw layer named 'contours14'. msContourLayerGenerateContour(): Unable to access file. GDALContourGenerate() failed: Input values and/or interval settings would lead to too many levels

[jmckenna]

with @rouault's changes the test case now passes, using @gsueur's coastline EXTENT values. Thanks Even!!

[gsueur]

I'm not sure it does, since I have seen a different behaviour from Mapserver (error message instead of crash) but not the expected one. Do you confirm the display of contours at high scale ? does the code to use is mapserver-7.6-branch ?
Best regards

[jmckenna]

@gsueur I'm testing using shp2img. Can you provide here new EXTENT values that cause your crash? Actually, I'll try to generate some...... (I used master branch)

[jmckenna]

Corsica island, zoomed way out, showing contours, with master branch.

[gsueur]

After comprehensive tests on branch master (and not 7.6 like the backport tag mislead me), I can confirm everything is fine.
Thanks for the help you both !