NetTopologySuite.IO.ShapeFile icon indicating copy to clipboard operation
NetTopologySuite.IO.ShapeFile copied to clipboard

Cannot read a polygon with a hole from a ShapeFile

Open dmitriytuma opened this issue 3 years ago • 2 comments

Hello! I have a problem with reading of one specific polygon in the attached shapefile (Example.zip) - this polygon has one hole, but it is read as MULTIPOLYGON EMPTY. Other polygons can be read without any problems.

But I can see the problem polygon, when I import the attached file in mapshaper, here it is marked in a red border: Problem polygon with a hole Meanwhile, there is another one polygon inside the hole

I have already tried types like ShapefileReader, ShapeDataReader, ShapefileDataReader, but result is the same.

One of my tries looks like:

using (var shapeFileDataReader = new ShapefileDataReader(shapeFilePathWithoutDotAndFileFormat, new GeometryFactory()))
{
    while (shapeFileDataReader.Read())
    {
        var geometry = shapeFileDataReader.Geometry;
    }
}

The projection of the shape file is UTM35N.

Please, can someone check this and help me?

Thank you!

dmitriytuma avatar Feb 26 '21 16:02 dmitriytuma

FInally I figured out what is a problem.

Both the problem polygon and its hole have IsCCW = true, which is really strange, because as I far as I know shells must be clockwise, and holes should be counterclockwise. But anyway - why and how mapshaper reads this incorrect polygon...

I would like to suggest a solution like this (in PolygonHandler.cs file, in Read method, after the initialization of lists shells and holes):

// Make the shell to be CW, and its holes - CCW.
if (shells.Count == 0 && holes.Any())
{
	holes.Sort(ProbeLinearRing);

	// Make the biggest "hole" as a shell.
	var newShell = holes.Last();
	shells.Add(newShell.IsCCW ? newShell.Reverse() as LinearRing : newShell);
	holes.RemoveAt(holes.Count - 1);

	// Validate the rest holes.
	for (int i = 0; i < holes.Count; i++)
	{
		holes[i] = holes[i].IsCCW ? holes[i] : holes[i].Reverse() as LinearRing;
	}
}

I am not sure how this will work in cases when we should have more than one shell - if someone can help me with this, it would be great! :)

Thank you!

dmitriytuma avatar Mar 02 '21 16:03 dmitriytuma

I tried to load the simble polygon "wrongly formatted" with arcmap successfully (see shell_bad_ccw.zip), using qgis I see that the polygon is visualized but there are problems with identify operation (i.e.: clicking on it). QGIS "WKT plugin" reads the data as MultiPolygon (((469380.33776138117536902 4860144.96628806926310062, 469413.05082988395588472 4860166.82644257601350546, 469433.65906083380104974 4860133.06535593792796135, 469405.23866025015013292 4860111.97170174960047007, 469380.33776138117536902 4860144.96628806926310062),(469884.03061507421080023 4859699.943781235255301, 469495.37857846909901127 4860354.11189078353345394, 469501.99340123386355117 4860329.19041431415826082, 469486.42985767225036398 4860285.03266813047230244, 469466.49032270145835355 4860272.67446962557733059, 469439.70877098350320011 4860259.814555699005723, 469420.69635292305611074 4860253.00041514821350574, 469400.66809175250818953 4860261.05747527442872524, 469359.1600068547995761 4860279.8071622671559453, 469324.76470285677351058 4860295.41577777080237865, 469303.8748967184801586 4860305.51915748789906502, 469256.13816838583443314 4860327.53448224626481533, 469237.84905035427073017 4860335.85752750653773546, 469203.60728499933611602 4860351.33781702537089586, 469190.26329871668713167 4860357.57371238246560097, 469189.8381122819846496 4860357.83147291745990515, 469138.53844545234460384 4860246.30646639782935381, 469118.43629730568500236 4860184.14369823317974806, 469118.48088543361518532 4860182.68971209414303303, 469122.25646422326099128 4860154.46180670894682407, 469126.35674297460354865 4860144.86082231905311346, 469126.92550657573156059 4860121.57886170502752066, 469223.67968696780735627 4860111.24572346825152636, 469298.08790045708883554 4860126.56679979432374239, 469304.02718955214368179 4860080.40053514018654823, 469546.97630507848225534 4859946.38098762370646, 469583.01080905337585136 4859928.52556280419230461, 469584.90992153104161844 4859927.3052488649263978, 469577.70299429236911237 4859864.34919862169772387, 469577.20574015128659084 4859862.31158902402967215, 469585.61254699126584455 4859860.33386100549250841, 469595.04299338022246957 4859858.83864029683172703, 469725.70676207466749474 4859838.04413766786456108, 469733.44323446520138532 4859773.57196296006441116, 469743.59146814869018272 4859689.32970975525677204, 469884.03061507421080023 4859699.943781235255301)))

Anyway, to me the input looks simply wrong: I will do a deeper check but the correct behavior to me is to throw an exception instead to read as MULTIPOLYGON EMPTY

EDIT: I just checked what happens inside PolygonHandler and there are A LOT of mess... basically:

  1. the shell is considered as a hole, because of the checks in the CCW
  2. the hole is considered as a shell, because of the same checks above
  3. the hole (actually, the shell interpreted as a hole) is removed and the returned result is a Polygon with a single shell and no holes POLYGON ((469380.33776138118 4860144.9662880693, 469413.05082988396 4860166.826442576, 469433.6590608338 4860133.0653559379, 469405.23866025015 4860111.97170175, 469380.33776138118 4860144.9662880693))

@dmitriytuma can you confirm the results? I don't see any MULTIPOLYGON EMPTY...

DGuidi avatar Jan 20 '22 14:01 DGuidi