Hello,

I decided to shift some of the geoprocessing functions from my vb6 application to ocx. They'll work faster there and may be usefull for community, especially for people with non-net languages.

Here is approximate list of what I plan to do:
- coordinates of lines/segments crossing;
- closest point and distance of shape to given point;
- centroid of polygon (from MapWinGeoProc);
- minimal distance between 2 shapes;
- shape rotation;

Relation between 2 shapes:
- intersection;
- boundries touch;
- common point;
- identical;
- within distance;
- inclusion.

Post here your thoughts on the topic.

Regards,
Sergei Leschinsky

Sergei,

Sounds great. Will you be making those functions from scratch or will you implement the GDAL/GEOS functions?

I'm very interested in the centroid of polygon function. The current one is calculating the center of gravity. If you have a L or C shaped polygon the centroid could be outside the polygon. Will you compensate for that?

I believe some of the functions you mention are already available in the MapWinGeoProc library. I do hope you are going to move them and redirect the MapWinGeoProc function call to the ocx and not copy them.

I think you've found it already but to be certain I think this would be a good location for your functions:
MapWinGIS:Utils

Thanks,

Paul


PS.
Your doing a great job. If you've got the time would you mind looking at this enhancement: [Labeling] PAL project

--
Don't forget to read the new documentation: www.mapwindow.org/documentation/mapwingis4.8
Join us Google+: MapWindow GIS Google+ Community
Join the MapWindow Group on LinkedIn! LinkedIn - MapWindow Group

Download the latest beta installer at:
tinyurl.com/mwMonthly 32-Bit
tinyurl.com/mwMonthlyx64 64-Bit
Follow me on Twitter MapWindow_nl to read when a new installer is published.

---
Paul Meems
The Netherlands
[www.bontepaarden.nl]
Release manager, configuration manager and
forum moderator of MapWindow GIS

Owner of MapWindow.nl - Support for
Dutch speaking users: www.mapwindow.nl

*******
Everything I say or write is my personal opinion and
not the opinion of the company I work for.
*******
View my profile on LinkedIn

Paul,

It seems you expect too much ;)

I just have vb6 functions which I wrote from scratch some months ago. They aren't very bright or flawless (my specialty is foresty after all). But they work already and it's not very difficult to shift them to cpp. Later on we may discuss and refine them.

I plan to add them to CUtils class and CShape class. I'll show API for discussion here.

Centroid is a copy from MapWinGeoProc (in single function though). I've already written and tested it. No compensation for L/C shaped polygons for now. But if you'll suggest how to do it (links/code or something) a can try.

What do you know about GDAL/GEOS functions? Can they be used to define relations between shapes (intersection,inclusion etc)?

Regards,
Sergei

About the centroid I couldn't get it to work either.
The PAL Library is using some algorithm to get the labels inside the polygon.
Perhaps it is wise to leave it for now and when we've added the PAL Library use its centroid function.

I don't know much about the GDAL/GEOS function but I do see a lot of interesting named files when I compile the library, like:
In SupportLibraries\GDAL\geos-cvs-x86\source\algorithm:
CentroidArea.cpp
LineIntersector.cpp
in SupportLibraries\GDAL\geos-cvs-x86\source\operation:
Buffer
Distance
Union
Valid

Useful links I found:
[gdal.org]
[n2.nabble.com]
[trac.osgeo.org]

Hope it helps,

Paul

--
Don't forget to read the new documentation: www.mapwindow.org/documentation/mapwingis4.8
Join us Google+: MapWindow GIS Google+ Community
Join the MapWindow Group on LinkedIn! LinkedIn - MapWindow Group

Download the latest beta installer at:
tinyurl.com/mwMonthly 32-Bit
tinyurl.com/mwMonthlyx64 64-Bit
Follow me on Twitter MapWindow_nl to read when a new installer is published.

---
Paul Meems
The Netherlands
[www.bontepaarden.nl]
Release manager, configuration manager and
forum moderator of MapWindow GIS

Owner of MapWindow.nl - Support for
Dutch speaking users: www.mapwindow.nl

*******
Everything I say or write is my personal opinion and
not the opinion of the company I work for.
*******
View my profile on LinkedIn

Hello,

I made an attemt to implement OGR/GEOS geoprocessing functions in ocx. I've completed functions for OGRGeometry / MapWinGis.Shape conversion.

IShape* CUtils::OGRGeometryToShape(OGRGeometry* oGeom)
OGRGeometry* CUtils::ShapeToOGRGeometry(IShape* shp)

I tested Shape-OGRGeometry-Shape conversion on shapefiles I have.

In all cases initial shape and converted shape have identical number of parts and points. Test on complete identity sometimes returns false for multipart polygons. For example for shapefile with complex polygons from bug 1068 ( [bugs.mapwindow.org] ) it returned false for 18/109 polys. But it's expected behaviour as OGR/GEOS can change the order of parts.

The time on two-way conversion is about 2 sec/3.5 polys (~20 points per poly).

Now I'm looking forward to introduce geoprocessing functions of OGR to ocx. Here is OGRGeometry class members: [gdal.org].

Your ideas on the topic are appreciated.

Regards,
Sergei

Hello,

Some questions emerged:

1.How can I test that geos_c.dll library is availible (in ocx)? I tested some new function and it gave me a crash first which disappeared after I moved geos_c.dll in Windows\System32 folder (the ocx is here also). Where can geos_c.dll be placed to be found in time?

2. Where geos_c.dll is used in ocx. It's said here [www.mapwindow.org] that in CUtils class but I don't see where for now.

Regards,
Sergei

P.S. I’ve implemented
STDMETHODIMP CUtils::ShapesRelation(IShape* Shape1, IShape* Shape2, tkSpatialRelation Relation, VARIANT_BOOL* retval)
and it’s working fine.

Hello,

1. Buffer function is ready and it's rather fast. It took 27 seconds to buffer about 3000 polygons shapes (~30 points each) without combining overlapping buffers. Distance 500 m, nQuadSegments = 30. It works equally good with other shape types.


2. I've encountered a problem while dealing with OGR/GEOS, which I'd like to document.

When debuging buffer function I got the next error:
HEAP[Forplan.exe]: Invalid Address specified to RtlFreeHeap( 014A0000, 0A929008 ) in immedeate window.

I've googled the following:

If you build a DLL or an EXE with debug linkage, it links to the debug version of new and delete. If you build one without debug, it links to the "release" versions. What may come as a surprise is that both versions can end up "successfully" linked into the same process space, as the EXE and each individual DLL get their own version!

[mail-archives.apache.org]
[social.msdn.microsoft.com]

Here is code which generate it.

GEOSGeom hGeosProduct = NULL;
// ...
if( hGeosProduct != NULL )
{
	size_t nSize = 0;
	unsigned char *pabyBuf = NULL;

	pabyBuf = GEOSGeomToWKB_buf( hGeosProduct, &nSize );	// memory is allocated inside GEOS
	if( pabyBuf == NULL || nSize == 0 ) return NULL;
	
	if( OGRGeometryFactory::createFromWkb( (unsigned char *) pabyBuf, NULL, &oRes, (int) nSize )!= OGRERR_NONE ) oRes = NULL;

	if( pabyBuf != NULL ) free( pabyBuf );		// here we free the memory (originally in ogrgeometryfactory.cpp
				                             // OGRGeometryFactory::createFromGEOS); error is here

	GEOSGeom_destroy( hGeosProduct );
	return oRes;
}

It's true - the problem doesn't appear with compiled dll. To solve it in debug mode it is needed:
- to free *pabyBuf inside geos_c.dll (there is no such possibility);
- maybe to use release version of gdal library (the free call is originally situated here)? do we use debug version of library now?

For now all calls which return new geometry from GEOS will produce this error in debug mode. For now it's:

BufferShape(Shape as MapWinGis.Shape, Distance as double, nQuadSegments as integer) as MapWinGis.Shape
ClipShape (Shape as MapWinGis.Shape, Clip as MapWinGis.Shape, oper as tkClipOperation) as MapWinGis.Shape
ConvexHull() as MapWinGis.Shape
Boundry() as MapWinGis.Shape

Any practical way to fix it? Chris, I'd like to hear your opinion.

Regards,
Sergei

Hi,

As my attempt with OGR/GEOS goes good, let's discuss API for new functins. First I'd like to introduce in ocx GEOS functions and maybe later my own functions which will be usefull still.

Below new fucntions that are in my list for now (I'll use vb6 notation for better understanding). I'd like to introduce all this functions to Shape class. I don't see any reason to shift them to Utils class as it was done with other functions before.

1)
Contains (Shape as MapWinGis.Shape) as boolean
Crosses (Shape as MapWinGis.Shape) as boolean
Disjoint (Shape as MapWinGis.Shape) as boolean
Equals (Shape as MapWinGis.Shape) as boolean
Intersects (Shape as MapWinGis.Shape) as boolean
Overlaps (Shape as MapWinGis.Shape) as boolean
Touches (Shape as MapWinGis.Shape) as boolean
Within (Shape as MapWinGis.Shape) as Boolean

This group can be unified with:
ShapesRelation (Shape1 as MapWinGis.Shape, Shape2 as MapWinGis.Shape, Relation as tkSpatialRelation) as boolean

2)
GetIntersection (Shape as MapWinGis.Shape) as MapWinGis.Shape
GetDifference (Shape as MapWinGis.Shape) as MapWinGis.Shape
GetUnion (Shape as MapWinGis.Shape) as MapWinGis.Shape
GetSymDifference (Shape as MapWinGis.Shape) as MapWinGis.Shape

This group can be unified with:
ClipShape (Shape as MapWinGis.Shape, Clip as MapWinGis.Shape, oper as tkClipOperation) as MapWinGis.Shape


Buffer (Shape as MapWinGis.Shape, Distance as double, nQuadSegments as integer) as MapWinGis.Shape
ConvexHull () as MapWinGis.Shape
Boundry () as MapWinGis.Shape

Centroid (Shape as MapWinGis.Shape) as MapWinGis.Point // MapWinGeoProc implementation would be faster I suppose
Distance (Shape1 as MapWinGis.Shape, Shape2 as MapWinGis.Shape) as double

Also I'd like to shift the following functions from Utils to Shape.
Length (Shape as MapWinGis.Shape) as double
Perimeter (Shape as MapWinGis.Shape) as double
Area (Shape as MapWinGis.Shape) as double
It's absolutely illogical to have these functions in Utils (Of course we must leave them there, but may have them in new location as well).

Description of OGRGeometry class members are here: [gdal.org]

So questions for discussion are:
1. Do you agree with naming?
2. Do you agree with location?
3. Do we need single functions, unified ones or both for first and second group?
4. Do we need Length, Perimeter and Area in Shape class as well as in Utils class?

Regards,
Sergei

Hi Sergei,

Naming and location looks good to me. It does make more sense having them on the Shape class, so I'd put them there. We do need to leave the old ones in Utils of course for backward compatibility, but otherwise looks great. Thanks for taking the time to do all of this.

Thanks!
--Chris

Sergei,

I agree with Chris. Your suggestions are fine.

About the debug/release issues. Since the last weekly build (yesterday) all supporting libraries are in release mode (including GDAL).
Hopefully that will solve your problem.

I saw you've added or will be adding a buffer method. Does that mean MW won't be using the NTS.Topology Buffer method anymore? That would be nice because it is giving trouble (bugs 1314, 436, 807 and 1349).

Thanks,

Paul

--
Don't forget to read the new documentation: www.mapwindow.org/documentation/mapwingis4.8
Join us Google+: MapWindow GIS Google+ Community
Join the MapWindow Group on LinkedIn! LinkedIn - MapWindow Group

Download the latest beta installer at:
tinyurl.com/mwMonthly 32-Bit
tinyurl.com/mwMonthlyx64 64-Bit
Follow me on Twitter MapWindow_nl to read when a new installer is published.

---
Paul Meems
The Netherlands
[www.bontepaarden.nl]
Release manager, configuration manager and
forum moderator of MapWindow GIS

Owner of MapWindow.nl - Support for
Dutch speaking users: www.mapwindow.nl

*******
Everything I say or write is my personal opinion and
not the opinion of the company I work for.
*******
View my profile on LinkedIn

Paul,

Thanks for info about supporting libraries.

About buffer:
I make changes to ocx only, but it'll be quite easy to redirect calls from MapWinGeoProc to the new function I think;

There are some issues we must bear in mind:
- testing is needed to catch possible bugs in new functions and to analyse cases when functions don't return resulting shape;
- GEOS Buffer function don't support the following enumerations used in MapWinGeoProc:
Enumerations.Buffer_HoleTreatment
Enumerations.Buffer_CapStyle
so changes in MapWindow GUI will be needed.

Regards,
Sergei



Edited 1 time(s). Last edit at 08/02/2009 01:32PM by Sergei.

Hello,

There is one more problem for discussion. Intersection function of GEOS returns Geometry collections which can hold objects of different types. For example:

intersection of 2 polys can return points, lines and polys simultaneously;
intersection of 2 lines can return both lines and points.

Such geometry collection can't be represented with a single shape. And since there are a number of shapes a problem arises how to return them to a calling program.

Possible ways are (which I see):
1. To use SAFEARRAY - I don't know will it be working with shapes; it's used for SelectShapes routine so far. If anybody have knowledge about it, please share it with me.

2. Return only shapes of the given type , for example:
Shape.Itersection(Shape as MapWinGis.Shape, returnType as MapWinGis.ShpFileType)
I don't like this as calculations must be perfomed several times to extract all results and time consumption is critical for this function.

3. To make 3 methods in utils class :
Utils.Intersection(Shape1, Shape2) as boolean // true if we have resulting shapes
Utils.GetNextResultShape() as MapWinGis.shape // return NULL whan shapes are over
Utils.ResultsCount() as long //
Result shapes are stored in utils class and are deleted on the next call of Intersection. This alternative I like more than others.

Please post, if you have ideas on the problem. Do you agree with API with 3 new functions in CUtils as shown above?

Problem with Invalid Address specified to RtlFreeHeap persists after I updated support libraries so any help is appreciated also.

Regards,
Sergei



Edited 1 time(s). Last edit at 08/05/2009 03:25AM by Sergei.

My preference would be to use a safearray - that way it looks straightforward to users, since they get an array of shapes back.

On the RtlFreeHeap thing, I have no idea I'm afraid.

Cheers,
--Chris

Hello,

I've commited changes.

There are following differences with the previous list:
1) Shape.IsValid() property (OGR function);
2) I didn't make individual methods for clip operations (difference, union, etc)
3) I added ShapesIntersection and NextResultShape methods to CUtils class. Their use (vb6):

  If Utils.ShapesIntersection(sf1.Shape(i), sf2.Shape(j)) Then
    Do
      Set shRes = Utils.NextResultShape
      If not shRes Is Nothing Then
        ' do something
      Else
        Exit Do
      End If
    Loop
  End If

Chris, maybe I'll rewrite them using SAFEARRAY. For now I've already implemented this approach.

The CShape.Clip, CUtils.ShapesIntersection and CUtils.NextResultShape are in test mode for now and their behaviour can be changed.

For now there is no check for geos_c.dll availibility. If there is no geos_c.dll you'll experience a crash most likely.

CShape.Centroid doesn't use external libraries.

Length, Perimeter, Area properties are shifted from CUtis. Call from Cutils are redirected to CShape class to avoid duplication of code.

I've tried to test most of the functions, but I expect there are still problems ahead. So let's test these bunch together.

Paul, add Ogr2MWshape.cpp and Ogr2MWshape.h to the project please (CShape and CUtils classes depend upon them now).

Regards,
Sergei

Sergei,

Could you add a brief list of all new features you've added and how to use them (pseudocode/VB6). I will then try to create some C# test scripts.
Will you also update the wiki, with your new functions so other can understand them as well and can use them. I'll add my C# code later as sample code.

Thanks,

Paul

--
Don't forget to read the new documentation: www.mapwindow.org/documentation/mapwingis4.8
Join us Google+: MapWindow GIS Google+ Community
Join the MapWindow Group on LinkedIn! LinkedIn - MapWindow Group

Download the latest beta installer at:
tinyurl.com/mwMonthly 32-Bit
tinyurl.com/mwMonthlyx64 64-Bit
Follow me on Twitter MapWindow_nl to read when a new installer is published.

---
Paul Meems
The Netherlands
[www.bontepaarden.nl]
Release manager, configuration manager and
forum moderator of MapWindow GIS

Owner of MapWindow.nl - Support for
Dutch speaking users: www.mapwindow.nl

*******
Everything I say or write is my personal opinion and
not the opinion of the company I work for.
*******
View my profile on LinkedIn

Paul,

the list of features was in this topic some posts ago. It differs a bit from final API, but I listed changes in previous post. The usage of functions is rather intuitive. Anyway you can use GDAL/OGR API reference for help [gdal.org].

I'll add description to wiki when I have time (it's rather difficult with it now).

Regards,
Sergei

Sergei,

I'm working on a script to show the use of the new GEOS functions
For testing purposes I create in code a shape.
When I create a point shape the IsValid returns True, but when creating a polygon it returns false.
That's ok, but it would be very nice if the Validation Error would be added to the LastError so I can see why it is not valid.
This is my code:

    private static MapWinGIS.Shape createTestShape(string scriptName, ShpfileType shpFileType)
    {
        // Create new shape:
        MapWinGIS.Shape testShape = null;

        try
        {
            ReportMessage(scriptName, "Make test shape");
            testShape = new MapWinGIS.Shape();
            testShape.Create(shpFileType);
            MapWinGIS.Point pt = new Point();
            // Centre of USA:
            pt.x = -109.5;
            pt.y = 37.5;
            int numPoint = testShape.numPoints;
            testShape.InsertPoint(pt, ref numPoint);
            if (shpFileType != ShpfileType.SHP_POINT)
            {
                // Add more points:
                pt = new Point();
                // Centre of USA:
                pt.x = -108.5;
                pt.y = 37.5;
                numPoint = testShape.numPoints;
                testShape.InsertPoint(pt, ref numPoint);

                // Add more points:
                pt = new Point();
                // Centre of USA:
                pt.x = -108.5;
                pt.y = 36.5;
                numPoint = testShape.numPoints;
                testShape.InsertPoint(pt, ref numPoint);

                // Add more points:
                pt = new Point();
                // Centre of USA:
                pt.x = -109.5;
                pt.y = 36.5;
                numPoint = testShape.numPoints;
                testShape.InsertPoint(pt, ref numPoint);

                if (shpFileType == ShpfileType.SHP_POLYGON)
                {
                    pt = new Point();
                    // Centre of USA:
                    pt.x = -109.5;
                    pt.y = 37.5;
                    numPoint = testShape.numPoints;
                    testShape.InsertPoint(pt, ref numPoint);
                }
            }
            ReportMessage(scriptName, "Validating the created shape: " + testShape.IsValid.ToString());
            ReportMessage(scriptName, "Validation error: " + testShape.get_ErrorMsg(testShape.LastErrorCode));
        }
        catch (Exception ex)
        {
            ReportMessage(scriptName, "Error createTestShape " + testShape.get_ErrorMsg(testShape.LastErrorCode) + " Exception: " + ex.ToString());
        }
        ReportMessage(scriptName, "NumPoints of testShape: " + testShape.numPoints.ToString());
        return testShape;
    }

ReportMessage(scriptName, "Validation error: " + testShape.get_ErrorMsg(testShape.LastErrorCode));

returns Validation error: No Error

Adding an invalid shape to the buffer function returns a null object. So I really need to know why it is invalid.

Thanks,

Paul

--
Don't forget to read the new documentation: www.mapwindow.org/documentation/mapwingis4.8
Join us Google+: MapWindow GIS Google+ Community
Join the MapWindow Group on LinkedIn! LinkedIn - MapWindow Group

Download the latest beta installer at:
tinyurl.com/mwMonthly 32-Bit
tinyurl.com/mwMonthlyx64 64-Bit
Follow me on Twitter MapWindow_nl to read when a new installer is published.

---
Paul Meems
The Netherlands
[www.bontepaarden.nl]
Release manager, configuration manager and
forum moderator of MapWindow GIS

Owner of MapWindow.nl - Support for
Dutch speaking users: www.mapwindow.nl

*******
Everything I say or write is my personal opinion and
not the opinion of the company I work for.
*******
View my profile on LinkedIn

Paul,

I didn't understand, if the function returns false for all polygons or the definite ones only? As long as I remember I tested poly shapefile and IsValid returned true for major part of the shapes.
If only some polygons are invalid, then send me those as a shapefile and I'll check what wrong with them.

Also I don't see that you added a part to the poly in your code. And it's perfect reason for not passing the conversion. Does your poly has a part?

Regards,
Sergei

Sergei,

My code is creating either a point, a line or a polygon.
The polygon is a square with no parts and 5 points, the last one equal to the first to close. Do I need a part? I haven't done that in the past.
The line is the same square but missing the last edge looking like a 'C'.
Tomorrow I'll start making a test script for the IsValid() method. To use some more data.

Any change of putting the validation error in the lastError of the shape?

Thanks,

Paul

--
Don't forget to read the new documentation: www.mapwindow.org/documentation/mapwingis4.8
Join us Google+: MapWindow GIS Google+ Community
Join the MapWindow Group on LinkedIn! LinkedIn - MapWindow Group

Download the latest beta installer at:
tinyurl.com/mwMonthly 32-Bit
tinyurl.com/mwMonthlyx64 64-Bit
Follow me on Twitter MapWindow_nl to read when a new installer is published.

---
Paul Meems
The Netherlands
[www.bontepaarden.nl]
Release manager, configuration manager and
forum moderator of MapWindow GIS

Owner of MapWindow.nl - Support for
Dutch speaking users: www.mapwindow.nl

*******
Everything I say or write is my personal opinion and
not the opinion of the company I work for.
*******
View my profile on LinkedIn

Paul,

Adding a part is a desired thing. Many functions in ocx process shape parts in first cycle and shape points in the nested cycle. If there is no part then a shape won't be processed. Ocx checks shapes and automatically adds one part for polyline and polygon after editing session (if I remember it right). Because of it adding shape without a part is working.

The best solution is to add a part for lines and polys in Shape.Create function. But it wasn't done, and we can't do it without harming existing code now.

I wrote strict check in GEOS functions, when a line or poly without a single part isn't processed further. It's possible to add single part automatically, but I don't like the idea. The shape won't be dispalyed right if there is more than one part and none of them were inserted. So I think it's necessary to return specific error code as you suggest.

Regards,
Sergei