Re: [postgis-users] Snapping linestrings with given tolerance

2012-11-20 Thread Guillaume Drolet
Hi Nicolas,

Thanks a lot for the tips. I will have a look at the thread and method you
suggested and keep you posted.

G
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


[postgis-users] Subject: Re: raster2pgsql 2.1.1 on Windows 7 64bit - ERREUR: erreur de syntaxe

2014-04-17 Thread Guillaume Drolet
Regina, Remi:

No space in filenames. GDAL version is 1.11dev. I don't have any
problems importing shapefiles into PostGIS.

Encoding was already UTF8 but I tried your suggestion (SET
PGCLIENT...) and I still have errors with SQL files having incomplete
last line.

I noticed that my LC_COLLATE and LC_TYPE are both set
'French.Canadian.1252'. Maybe this is the source of the problem. I
tried creating a new cluster with --no-locale to make some tests with
raster2pgsql but I screwed up and messed with my current cluster in
the process...

Will do more tests when I'm back from Easter Holiday and I fixed my
cluster problem.
Thanks for your help.

Best,

Guillaume
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


[postgis-users] Subject: Re: raster2pgsql 2.1.1 on Windows 7 64bit - ERREUR: erreur de syntaxe

2014-04-24 Thread Guillaume Drolet
I tried the same command as in my previous post on my PC at home. It's
Windows 7 64bit as well and the cluster has the same locale
('French.Canadian.1252'.) and encoding (UTF8).

The command was successful so my work colleague and I, who experiences
the same behaviour as the one described previously on this list,
suspect that it may be due to some problems or conflict between
PostGIS and the hardware (we have similar systems): is this a sensible
hypothesis?

Thanks for your help,

Guillaume

> Message: 4
> Date: Thu, 17 Apr 2014 17:03:58 -0400
> From: Guillaume Drolet 
> To: "postgis-users@lists.osgeo.org" 
> Subject: [postgis-users] Subject: Re: raster2pgsql 2.1.1 on Windows 7
> 64bit - ERREUR: erreur de syntaxe
> Message-ID:
> 
> Content-Type: text/plain; charset=UTF-8
>
> Regina, Remi:
>
> No space in filenames. GDAL version is 1.11dev. I don't have any
> problems importing shapefiles into PostGIS.
>
> Encoding was already UTF8 but I tried your suggestion (SET
> PGCLIENT...) and I still have errors with SQL files having incomplete
> last line.
>
> I noticed that my LC_COLLATE and LC_TYPE are both set
> 'French.Canadian.1252'. Maybe this is the source of the problem. I
> tried creating a new cluster with --no-locale to make some tests with
> raster2pgsql but I screwed up and messed with my current cluster in
> the process...
>
> Will do more tests when I'm back from Easter Holiday and I fixed my
> cluster problem.
> Thanks for your help.
>
> Best,
>
> Guillaume
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


[postgis-users] Subject: Re: raster2pgsql 2.1.1 on Windows 7 64bit - ERREUR: erreur de syntaxe

2014-04-28 Thread Guillaume Drolet
Hi Mateusz,

Thanks for helping. What you suggested, i.e. outputting raster2pgsql
to an SQL file instead of piping it to psql, I have done already. I
also posted examples of where the lines are broken in previous post
(see also Jean-Daniel Sylvain's older post on the same topic).

It doesn't seem to break on non-ASCII characters. For example, this
one breaks in the middle of the file name at line 168:

...
INSERT INTO "mod09a1.a249.h12v04.005.2006268184041_sur_refl_b01"
("rast","filename") VALUES
('010100A544B7031A11713FA544B7031A1171BF86BE196C857752C095AFEE8CEE6E4740E61032001F0087004D4F44303941312E41323030303034392E6831327630342E3030352E323030363236383138343034315F7375725F7265666C5F6230312E74696600'::raster,'MOD09A1.A249.h12v04.005.2006268184041_sur_refl_b01.tif');
INSERT INTO "mod09a1.a249.h12v04.005.2006268184041_sur_refl_b01"
("rast","filename") VALUES
('010100A544B7031A11713FA544B7031A1171BF59D7C60F306A52C095AFEE8CEE6E4740E61032001F0087004D4F44303941312E41323030303034392E6831327630342E3030352E323030363236383138343034315F7375725F7265666C5F6230312E74696600'::raster,'MOD09A1.A249.h12v04.005.2006268184041_sur_refl_b01.tif');
INSERT INTO "mod09a1.a249.h12v04.005.2006268184041_sur_refl_b01"
("rast","filename") VALUES
('010100A544B7031A11713FA544B7031A1171BF2BF073B3DA5C52C095AFEE8CEE6E4740E61032001F0087004D4F44303941312E41323030303034392E6831327630342E3030352E323030363236383138343034315F7375725F7265666C5F6230312E74696600'::raster,'MOD09A1.A249.h12v04.005.2006268184041_sur_refl_b01.tif');
INSERT INTO "mod09a1.a249.h12v04.005.20062

Or, this one (same tif file; I only shortened the filename) breaks in
the middle of a cell (or tile?) representation, still at line 168:

...
INTO "mod09a1.a249.sur_refl_b01" ("rast","filename") VALUES
('010100A544B7031A11713FA544B7031A1171BF86BE196C857752C095AFEE8CEE6E4740E61032001F0087004D4F44303941312E41323030303034392E7375725F7265666C5F6230312E74696600'::raster,'MOD09A1.A249.sur_refl_b01.tif');
INSERT INTO "mod09a1.a249.sur_refl_b01" ("rast","filename") VALUES
('010100A544B7031A11713FA544B7031A1171BF59D7C60F306A52C095AFEE8CEE6E4740E61032001F0087004D4F44303941312E41323030303034392E7375725F7265666C5F6230312E74696600'::raster,'MOD09A1.A249.sur_refl_b01.tif');
INSERT INTO "mod09a1.a249.sur_refl_b01" ("rast","filename") VALUES
('010100A544B7031A11713FA544B7031A1171BF2BF073B3DA5C52C095AFEE8CEE6E4740E61032001F0087004D4F44303941312E41323030303034392E7375725F7265666C5F6230312E74696600'::raster,'MOD09A1.A249.sur_refl_b01.tif');
INSERT INTO "mod09a1.a249.sur_refl_b01" ("rast","filename") VALUES
('010100A544B7031A11713FA544B7031A117

This is VERY frustrating and I hope an answer will be found!

Guillaume


Message: 2
Date: Thu, 24 Apr 2014 21:43:43 +0200
From: Mateusz ?oskot 
To: PostGIS Users Discussion 
Subject: Re: [postgis-users] Subject: Re: raster2pgsql 2.1.1 on
Windows 7 64bit - ERREUR: erreur de syntaxe
Message-ID:

Content-Type: text/plain; charset=UTF-8

On 24 April 2014 21:19, Guillaume Drolet  wrote:
> I tried the same command as in my previous post on my PC at home. It's
> Windows 7 64bit as well and the cluster has the same locale
> ('French.Canadian.1252'.) and encoding (UTF8).
>
> The command was successful so my work colleague and I, who experiences
> the same behaviour as the one described previously on this list,
> suspect that it may be due to some problems or conflict between
> PostGIS and the hardware (we have similar systems): is this a sensible
> hypothesis?

I doubt it. I use raster2pgsql frequently on Windows, since PostGIS pre-2.0
and I have never experienced any such issues.

raster2pgsql is a standalone command line utility that generates
SQL commands and outputs them to stdout.
raster2pgsql has no run-time dependency on PostgreSQL or PostGIS.

First thing to check is: does raster2pgsql generate valid output,
if redirected to file, does the SQL commands in output file
look well-formed and complete?

If not, which SQL lines are broken?
Do the broken lines correspond to any character data in non-ASCII encodings?
For example, check if you run raster2pgsql with -F switch to add
column with the name of
your raster file and the raster file name contains characters with
French accents, etc.

I'd suggest, don't use psql.exe or load any data to PostgreSQL until
you confirm raster2pgsql generates text output with complete well-formed and
valid SQL commands first.

Best regards,
--
Mateusz  ?oskot, http://mateusz.loskot.net
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


[postgis-users] Subject: Re: raster2pgsql 2.1.1 on Windows 7 64bit - ERREUR: erreur de syntaxe

2014-04-30 Thread Guillaume Drolet
Hi Regina and others:

We found the problem with raster2pgsql! It was my version of GDAL
(1.11dev) that my GDAL environment variable was pointing to that
created the bug. I simply moved my GDAL folder to some other location
(without changing GDAL environment variable) and ran raster2pgsql
again so that PostGIS could use its own GDAL instead. This resulted in
a complete (i.e. no broken line at the end) output SQL file from
raster2pgsql.

Since I want to keep using my GDAL version (1.11dev, with all my
add-ons compiled and the version used by other software), this raises
the question: how to force PostGIS to not use the GDAL installation
referenced by the GDAL environment variable but instead use its own
GDAL?

Thanks!

G

2014-04-30 15:00 GMT-04:00  :
> Send postgis-users mailing list submissions to
> postgis-users@lists.osgeo.org
>
> To subscribe or unsubscribe via the World Wide Web, visit
> http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users
> or, via email, send a message with subject or body 'help' to
> postgis-users-requ...@lists.osgeo.org
>
> You can reach the person managing the list at
> postgis-users-ow...@lists.osgeo.org
>
> When replying, please edit your Subject line so it is more specific
> than "Re: Contents of postgis-users digest..."
>
>
> Today's Topics:
>
>1. Re: Subject: Re: raster2pgsql 2.1.1 on Windows 7  64bit -
>   ERREUR: erreur de syntaxe (Paragon Corporation)
>
>
> --
>
> Message: 1
> Date: Tue, 29 Apr 2014 18:51:16 -0400
> From: "Paragon Corporation" 
> To: "'PostGIS Users Discussion'" 
> Subject: Re: [postgis-users] Subject: Re: raster2pgsql 2.1.1 on
> Windows 7   64bit - ERREUR: erreur de syntaxe
> Message-ID: <40C3ECB89B9F4C388A6CE0F390D877C5@O>
> Content-Type: text/plain;   charset="us-ascii"
>
> Guillaume,
>
> Did you say this works fine on your home pc, but not your work one?
>
> I'm wondering if maybe its some malware or virus checker shutting down the
> execution of raster2pgsql for some reason.
>
> Have you looked at the windows event viewer to see if it is showing any
> messages?  I recall one time a behind the scenes VMWare virus checker kept
> on deleting some of my Msys files because it thought it was hacker ware. No
> duh. Had to tell sysadmins to exclude the folder from their checking.
>
> Hope that helps,
> Regina
>
> -Original Message-
> From: postgis-users-boun...@lists.osgeo.org
> [mailto:postgis-users-boun...@lists.osgeo.org] On Behalf Of Guillaume Drolet
> Sent: Monday, April 28, 2014 2:30 PM
> To: postgis-users@lists.osgeo.org
> Subject: [postgis-users] Subject: Re: raster2pgsql 2.1.1 on Windows 7 64bit
> - ERREUR: erreur de syntaxe
>
> Hi Mateusz,
>
> Thanks for helping. What you suggested, i.e. outputting raster2pgsql to an
> SQL file instead of piping it to psql, I have done already. I also posted
> examples of where the lines are broken in previous post (see also
> Jean-Daniel Sylvain's older post on the same topic).
>
> It doesn't seem to break on non-ASCII characters. For example, this one
> breaks in the middle of the file name at line 168:
>
> ...
> INSERT INTO "mod09a1.a249.h12v04.005.2006268184041_sur_refl_b01"
> ("rast","filename") VALUES
> ('010100A544B7031A11713FA544B7031A1171BF86BE196C857752C095AFEE8CEE6E4740
> E61032001F0087004D4F44303941312E
> 41323030303034392E6831327630342E3030352E323030363236383138343034315F7375725F
> 7265666C5F6230312E74696600'::raster,'MOD09A1.A249.h12v04.005.20062681840
> 41_sur_refl_b01.tif');
> INSERT INTO "mod09a1.a249.h12v04.005.2006268184041_sur_refl_b01"
> ("rast","filename") VALUES
> ('010100A544B7031A11713FA544B7031A1171BF59D7C60F306A52C095AFEE8CEE6E4740
> E61032001F0087004D4F44303941312E
> 41323030303034392E6831327630342E3030352E323030363236383138343034315F7375725F
> 7265666C5F6230312E74696600'::raster,'MOD09A1.A249.h12v04.005.20062681840
> 41_sur_refl_b01.tif');
> INSERT INTO "mod09a1.a249.h12v04.005.2006268184041_sur_refl_b01"
> ("rast","filename") VALUES
> ('010100A544B7031A11713FA544B7031A1171BF2BF073B3DA5C52C095AFEE8CEE6E4740
> E61032001F0087004D4F44303941312E
> 41323030303034392E6831327630342E3030352E323030363236383138343034315F7375725F
> 7265666C5F6230312E74696600'::raster,'

[postgis-users] raster2pgsql 2.1.1 on Windows 7 64bit - ERREUR: erreur de syntaxe

2014-05-02 Thread Guillaume Drolet
Regina,

When I do

SELECT postgis_full_version();

It says (from what I remember. I won't have access to my database
until next week) :

GDAL 1.11dev (of course if my 1.11dev installation is at the place
described by the GDAL variable). If my GDAL installation is moved as I
did in my previous post, the GDAL version shown by SELECT
postgis_full_version(); is 1.10.

So it looks like a problem with raster2pgsql and rather than with
PostGIS as the loading of the SQL file from raster2pgsql runs fine. I
don't know who maintains this program but maybe it'd be worth checking
it against versions of GDAL after 1.10.

Regards,

Guillaume

>
> Message: 2
> Date: Wed, 30 Apr 2014 21:24:06 -0400
> From: "Paragon Corporation" 
> To: "'PostGIS Users Discussion'" 
> Subject: Re: [postgis-users] Subject: Re: raster2pgsql 2.1.1 on
> Windows 7   64bit - ERREUR: erreur de syntaxe
> Message-ID: 
> Content-Type: text/plain;   charset="us-ascii"
>
> Guillaume,
>
> It was my understanding that PostGIS only uses the libgdal in the same
> folder as raster2pgsql since that would be the first in its path. However it
> uses the environment variable for the proj files and so forth.  So I'm still
> a bit surprised this is an issue unless those files are vastly different
> from 1.10.
>
> When you do
>
> SELECT postgis_full_version();
>
> What does it give for the GDAL version.
>
>
> If its just a during load problem, I think you can just simply override the
> environment variable in your loader batch script
>
> SET GDAL_DATA="C:\Program Files\PostgreSQL\9.3\gdal-data"
>
> Of course this won't fix your running install.  In PostGIS 2.2 we introduced
> a PostgreSQL GUC variable to handle this:
>
> http://postgis.net/docs/manual-dev/postgis_gdal_datapath.html
>
> But unfortunately that won't help you for 2.1 and not sure how much that
> helps with raster2pgsql anyway since I don't think raster2pgsql makes
> connection to Postgres to be able to even see that variable.
>
> What I could do is create a build that uses the newer GDAL to see if that
> improves your situation (you can be my test guinea pig :) ).  IN 2.2 I plan
> to ship with GDAL 1.11 and probably upgrade next  2.1 to GDAL 1.11 as well
> since unfortunately they have to share the same GDAL on same PostgreSQL
> install and some people I expect will be running both.
>
>
> Opinions?
> Regina
> http://www.postgis.us
>
>
> -Original Message-
> From: postgis-users-boun...@lists.osgeo.org
> [mailto:postgis-users-boun...@lists.osgeo.org] On Behalf Of Guillaume Drolet
> Sent: Wednesday, April 30, 2014 4:39 PM
> To: postgis-users@lists.osgeo.org
> Subject: [postgis-users] Subject: Re: raster2pgsql 2.1.1 on Windows 7 64bit
> - ERREUR: erreur de syntaxe
>
> Hi Regina and others:
>
> We found the problem with raster2pgsql! It was my version of GDAL
> (1.11dev) that my GDAL environment variable was pointing to that created the
> bug. I simply moved my GDAL folder to some other location (without changing
> GDAL environment variable) and ran raster2pgsql again so that PostGIS could
> use its own GDAL instead. This resulted in a complete (i.e. no broken line
> at the end) output SQL file from raster2pgsql.
>
> Since I want to keep using my GDAL version (1.11dev, with all my add-ons
> compiled and the version used by other software), this raises the question:
> how to force PostGIS to not use the GDAL installation referenced by the GDAL
> environment variable but instead use its own GDAL?
>
> Thanks!
>
> G
>
> 2014-04-30 15:00 GMT-04:00  :
>> Send postgis-users mailing list submissions to
>> postgis-users@lists.osgeo.org
>>
>> To subscribe or unsubscribe via the World Wide Web, visit
>> http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users
>> or, via email, send a message with subject or body 'help' to
>> postgis-users-requ...@lists.osgeo.org
>>
>> You can reach the person managing the list at
>> postgis-users-ow...@lists.osgeo.org
>>
>> When replying, please edit your Subject line so it is more specific
>> than "Re: Contents of postgis-users digest..."
>>
>>
>> Today's Topics:
>>
>>1. Re: Subject: Re: raster2pgsql 2.1.1 on Windows 7  64bit -
>>   ERREUR: erreur de syntaxe (Paragon Corporation)
>>
>>
>> --
>>
>> Message: 1
>> Date: Tue, 29 Apr 2014 18:51:16 -0400
>> From: "Paragon Corporation" 
>> To: "'PostGIS Users Discussion'&quo

[postgis-users] Topology: cannot delete slivers (or gaps)

2014-10-31 Thread Guillaume Drolet
Hi!

I'm am turning to you all hoping you'll have some tricks to help with
this problem:

I created a topology for 2290 polygons and populated a topogeom column
using these commands:

SELECT topology.CreateTopology('de_20k_topo', find_srid('syshiera',
'de_20k', 'geom_32198'), 1e-4);
SELECT topology.AddTopoGeometryColumn('de_20k_topo', 'syshiera',
'de_20k', 'topogeom', 'MULTIPOLYGON');

DO $$DECLARE r record;
BEGIN
  FOR r IN SELECT * FROM syshiera.de_20k LOOP
BEGIN
  UPDATE syshiera.de_20k SET topogeom = toTopoGeom(geom_32198,
'de_20k_topo', 1, 1e-4)
  WHERE gid = r.gid;
EXCEPTION
  WHEN OTHERS THEN
RAISE WARNING 'Loading of record % failed: %', r.gid, SQLERRM;
END;
  END LOOP;
END$$;

-- For testing:
SELECT topology.CopyTopology('de_20k_topo', 'de_20k_topo_tests');


I need to fix 1) extra nodes along edges that are not at an
intersection and 2) sliver polygons.

For cases in 1, I was successful using 'ST_NewEdgeHeal' and a plPgsql script.

For cases in 2, I can remove some slivers with 'ST_RemEdgeModFace' but
some of them I just can't delete:

I get either:

ERROR: TopoGeom 1123 in layer 1 (de_20k_topo_tests.LAYER1.) cannot be
represented healing faces 1514 and 1502

Or, when trying with remove a node with 'topology.ST_NewEdgeHeal':

ERROR:  SQL/MM Spatial exception - other edges connected (3912)

In this particular example, two edges are connected two the same two
nodes, similar to what was discussed here:

http://postgis.17.x6.nabble.com/PostGIS-1955-Exception-when-2-edges-passed-to-ST-ModEdgeHeal-that-are-attached-to-the-same-2-nodes-td4999383.html

But according to that discussion thread, this was fixed in an earlier
version (I'm using topology version 2.1.3).

Maybe my problem is different but I'll need your help to get it
sorted. Here's an example showing the small face (1502) I want remove:

https://dl.dropboxusercontent.com/u/5196336/example1.bmp

Thanks for your help!
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


Re: [postgis-users] Topology: cannot delete slivers (or gaps)

2014-10-31 Thread Guillaume Drolet
Hi Sandro,

Thanks for your help!

topology.ValidateTopology gives this:

error  id1   id2
character varying integer integer
--   ---   
face has no rings 1831

If I calculate the area of face 1831
(st_area(st_getfacegeometry('de_20k_topo', face_id))), I get a NULL value.

I also find I have some 650 of these slivers with an area below 75 m2 while
most of my correct polygons have an area above 300 000 000 m2.

Re the edges passed to ST_NewEdgeHeal, in the example I sent it was edges
3908 and 3909.

If you think about something else that could be problematic, let me know.
Meanwhile I'll dig into the approach you suggested in your reply (i.e.
fixing TopoGeometry objects).

Thanks!

PS. Forgot to mention that in the first topology I had created, I used a
tolerance of 2.0 meters instead of 1e-4 meters. While that yielded much less
slivers, I was getting way more of those nodes not located at an
intersection. What is best?





--
View this message in context: 
http://postgis.17.x6.nabble.com/Topology-cannot-delete-slivers-or-gaps-tp5007250p5007257.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


Re: [postgis-users] Topology: cannot delete slivers (or gaps)

2014-11-03 Thread Guillaume Drolet
Sandro Santilli wrote
>> I suggest you take a look at edges having face 1831 on either
>> the right or left edge. Do any edge exist ? 
> 
> edge 5187 has face 1831 as its right face
> 
>> The correct procedure to drop those slivers would be to assign
>> those small areas to the (closest|biggest|leftmost|youpickit)
>> TopoGeometry object, to one and only one, and then, if you really
>> need it, drop the edge internal to the TopoGeometry.
> 
>> I made a QGIS plugin which may help you do some of these things
>> manually, which is useful in some cases, mainly to understand the
>> issues.
> 
> The plugin you made is a great tool (thanks). I will first try to get the
> number of slivers down to 
> a manageable number, which I had with setting the tolerance at 2.0 m, and
> then apply the
> TopoGeometry approach above.
> 
>> So the error was right about there being another edge attached to the
>> node:
> 
>>  ERROR:  SQL/MM Spatial exception - other edges connected (3912)
> 
>> I'm guessing 3912 and 3908 are the two parallel edges in the picture
>> you attached: https://dl.dropboxusercontent.com/u/5196336/example1.bmp
> 
> You're right, they are the two parallel edges. And it's the same for most
> of
> the some 650 slivers I have.
> 
>> Again, QGIS may help you getting more info out of the visual, enable
>> the DBManager core plugin, select your topology schema and select
>> Schema->TopoViewer from the menu.
> 
> I didn't know about the TopoViewer tool. I was importing the topology
> tables (node, edge, face)
> like any other PG geometry tables. Will explore with this new tool.
> 
> 
>> It's probably easier to clean up nodes for an areal topology than to
>> clean up faces. With my QGIS plugin (PostGIS Topology Editor) you can
>> select all nodes and click on a button which will drop all the ones
>> that are safe to drop.
> 
> You're right! I'll create my topology again but with the 2.0 m tolerance,
> then run
> my function for getting rid of nodes that are not at intersections and
> finally,
> I will use the approach you suggested of working with the TopoGeometry
> objects and the
> relation table to eliminate the remaining slivers.
> 
> Will keep you posted on how it goes. Again, I really appreciate your help.
> Thanks.
> 
> Guillaume.





--
View this message in context: 
http://postgis.17.x6.nabble.com/Topology-cannot-delete-slivers-or-gaps-tp5007250p5007266.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


Re: [postgis-users] Topology: cannot delete slivers (or gaps)

2014-11-04 Thread Guillaume Drolet
Hi Sandro and others,

First, thanks Jon for the hint on pprepair: I ran pprepair on my shapefile
and I'm now in the process of 
importing it into PostGIs for topology creation. Then I'll be able to tell
if pprepair really fixed everything.

But because I'm stubborn and never give up, I still want to succeed on
fixing TopoGeometry objects in PostGIS (bear with me, what follows is a bit
long):

Sandro, you said:

> The fix here would be to find faces that take part in the composition
> of multiple TopoGeometry objects and those NOT taking part in the
> composition of any TopoGeometry object. I've seen GRASS GUIs referring
> to those 2 classes as "polygon0" and "polygon2", enough that I've been
> thinking it could be useful to have them shown as layers in the QGIS
> viewer for PostGIS Topology. You should be able to extract those faces
> with a query to the relation and the layers tables.

> So basically you first have to fix the TopoGeometry objects and then
> can safely drop the no-more-needed edges. 

First, I looked at faces NOT taking part in the composition of any
TopoGeometry object:

/SELECT face_id FROM de_20k_topo.face WHERE face_id NOT IN (SELECT
element_id

   
FROM de_20k_topo.relation);/
face_id
0
2418
2650
2663
2766
3077
3118
3146
3267
3270

None of these 10 faces are involved with my three remaining slivers... And
when I look at them, they ARE related to edges, as a next_face or a
right_face. Is there any useful information in that for my purpose?

Second, I tried to work along the lines suggested by Sandro, using a new
topology (tolerance 2.0 m). As in one of my previous attempts, I
successfully fixed 20 of the 23 slivers using ST_RemEdgeModFace. For the
remaining three slivers, however, I tried modyfying TopoGeometry objects but
with no success. I must admit I'm a bit lost and not too sure if I really
understand what TopoGeometry objects are:

For example, let's look at problematic face_id: 2227
(https://dl.dropboxusercontent.com/u/5196336/example_2227.bmp). In this
case, the light grey line is where a single edge should be. Blue edges 6381
and 6341 should not be there. The light grey line should start from the node
at the bottom right and connect at its other end somewhere in the middle of
edge 6340:

I thought, let's first remove edge 6346, which should not exist in the first
place, and then split edge 6340 with a node and connect it to node 4753 on
the lower right corner of the image:

SELECT ST_RemEdgeModFace('de_20k_topo', 6346); 

As in my previous post, this yields:

ERROR:  TopoGeom 1546 in layer 1 (syshiera.de_20k.topogeom) cannot be
represented healing faces 2225 and 2227
CONTEXT:  instruction SQL « SELECT topology._ST_RemEdgeCheck(toponame,
topoid, e1id, e1rec.left_face, e1rec.right_face) »
function PL/pgsql st_remedgemodface(character varying,integer), line 43 at
PERFORM

So I looked at TopoGeom 1546:

/SELECT * FROM de_20k_topo.relation WHERE topogeo_id = 1546;/

topogeo_id;layer_id;element_id;element_type
1546;1;2225;3
1546;1;2233;3

Edge 2227 is missing from TopoGeom 1546: Where is it?

/SELECT * FROM de_20k_topo.relation WHERE element_id = 2227;/

topogeo_id;layer_id;element_id;element_type
1551;1;2227;3

There it is, in TopoGeom 1551. What else is there:

/SELECT * FROM de_20k_topo.relation WHERE topogeo_id = 1551;/

topogeo_id;layer_id;element_id;element_type
1551;1;2225;3
1551;1;2234;3
1551;1;2227;3

So edges 2225 and 2227 are associated with TopoGeom 1551. Will change that:

/
UPDATE de_20k_topo.relation SET topogeo_id = 1546
WHERE topogeo_id = 1551 AND element_id = 2227;/

Query returned successfully: one row affected, 481 ms execution time.

And then again:

/SELECT ST_RemEdgeModFace('de_20k_topo', 6346); /

ERROR:  TopoGeom 1551 in layer 1 (syshiera.de_20k.topogeom) cannot be
represented healing faces 2225 and 2227
CONTEXT:  instruction SQL « SELECT topology._ST_RemEdgeCheck(toponame,
topoid, e1id, e1rec.left_face, e1rec.right_face) »
function PL/pgsql st_remedgemodface(character varying,integer), line 43 at
PERFORM

Obviously, as mentionned above, I don't quite understand how TopoGeometry
objects work and the link between TopoGeoms, geometries and primitives. I
need a bit of education here.

Thanks a lot for your patience and your help.

PS. Remi, thanks for suggesting GRASS but I've never used it and the
learning curve seems pretty steep, no?




--
View this message in context: 
http://postgis.17.x6.nabble.com/Topology-cannot-delete-slivers-or-gaps-tp5007250p5007273.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users

Re: [postgis-users] Topology: cannot delete slivers (or gaps)

2014-11-05 Thread Guillaume Drolet
Sandro Santilli wrote
>> Except face 0, which is the "universe" face, all the others
>> should be "gaps" in your topology. You can zoom to the face's
>> bounding box ("mbr" field) to take a closer look (do you have qgis
>> set up by now?).
> 
> Ok. I had already looked at some of them but as face-derived geometries
> (using ST_GetFaceGeometry) and not as real faces. As geometries, to me
> they look like polygons and I didn't realize they were in fact holes. See
> this example in QGIS for face_id 2418 turned into a polygon:
> https://dl.dropboxusercontent.com/u/5196336/geomFromFace_id_2418.bmp. Now,
> look at its mbr representation:
> https://dl.dropboxusercontent.com/u/5196336/mbrFace_id_2418.bmp. How do
> you tell it's a hole and not just a regular face?
> 
>> What you could do programmatically is for each face to fetch the list
>> of edges bounding it and pick the face on the other side as a candidate;
>> then for each candidate face you find TopoGeometry objects containing
>> it in their definition and pick one to assign that orphaned face
>> (you could use smaller TopoGeometry, or larger, or whatever).
> 
> Let's try this or face_id 2418 (which is, btw, a meaningful unit and not
> an error...):
/
> SELECT f.face_id, e.edge_id, e.left_face, e.right_face 
> FROM de_20k_topo.face f JOIN 
>  de_20k_topo.edge_data e ON (f.face_id = e.left_face OR
>f.face_id = e.right_face)
> WHERE face_id = 2418;
/
> 
*
> face_id, edge_id, left_face, right_face
*
> 2418;2663;2418;950
> 2418;6986;2418;340
> 2418;2672;2418;951
> 2418;2665;2418;344
> 2418;6989;2418;2417
> 
> I pick the first edge (2663) and take it's right_face, 950, as a
> candidate:
/
> SELECT topogeo_id 
> FROM de_20k_topo.relation
> WHERE element_id = 950;
/
> 
*
> topogeo_id
*
> 814
> 
> Only TopoGeom 814 has face_id 950 as one of its elements. Will assign
> face_id 2418
> to this TopoGeom:
/
> INSERT INTO de_20k_topo.relation 
> VALUES (814, 1, 2418, 3);
/
> 
> Then what? Is that it? Does that mean that that big hole 2418 is now
> unioned with face_id 950 (in which case it isn't what I want as it's a
> meaningful ecological unit)?
> 
> 
>> It might be better first to deal with the faces that are associated
>> with multiple TopoGeometry, just to avoid assigning an orphaned face
>> to a TopoGeometry based on it's containing a face which will later
>> be removed from the definition because duplicated.
> 
> Good point.
> 
>> Eh, that part is a still a bit weak in documantation, maybe.
>> A TopoGeometry is a Geometry defined by its components.
> 
>> In your case (areal TopoGeometry objects) it is a (multi)Polygon
>> defined by a set of faces. Note that QGIS is also able to
>> show you the TopoGeometry objects, and lets you edit them "spatially"
>> (less robust than editing them by adding/removing components but with
>> some care easier to use).
> 
>>> In this case, the light grey line is where a single edge should be. Blue
>>> edges 6381
>>> and 6341 should not be there. The light grey line should start from the
>>> node
>>> at the bottom right and connect at its other end somewhere in the middle
>>> of
>>> edge 6340:
> 
>> Uhm, edge identifiers are missing from that dump.
>> I take it 6381 is the lower blue and 6341 the upper blue.
> 
> You're correct here.
> 
>> I suggest you load the TopoGeometry layer with qgis and see for yourself.
> 
> When you say load the TopoGeometry layer, I assume you mean the topogeom
> column I populated with:
/
> UPDATE syshiera.de_20k SET topogeom = toTopoGeom(geom_32198,
> 'de_20k_topo', 1, 2.0)
/
> 
> Right?
> 
> When I load it in QGIS using the PostGIS table loader, I can visualize it
> and see there are holes. E.g.
> https://dl.dropboxusercontent.com/u/5196336/holeInTopogeometry.bmp
> 
>> Unfortunately (that'd be useful to add in qgis) it won't be easy to
>> find TopoGeometry 1546 as extracting the id of a TopoGeometry requires
>> running an expression on the database: id(topogeom).
> 
> I can either
/
> SELECT * 
> FROM syshiera.de_20k
> WHERE id(topogeom) = 1546;
/
> 
> or, as I mentioned before, I can superimpose my 'faces' layer, which I
> derived from de_20k_topo.face
> using ST_GetFaceGeometry, to 'topogeom' to get the topogeo_id. No?
> 
>> You could create a view selecting just that topoGeometry, or a subset
>> of them, or you could just load them all, hoping it won't be too
>> expensive
>> (another spot where qgis support for postgis topology could be improved).
> 
> I agree, just loading the TopoViewer freezes QGIS for a while and I don't
> lack ressources on my machine.
> 
>> Chances are it would be enough to assign face 2227 to the TopoGeometry
>> with id 1546 (likely the one not-assigned). But as face 2227 was not
>> in your initial set (unassigned faces) I guess there's another
>> TopoGeometry
>> which does have 2227 in his definition, but not the one on the right.
>> In that case you have to pick which of the two contending TopoGeometry
>> objects should get face 2227 assign

Re: [postgis-users] Topology: cannot delete slivers (or gaps)

2014-11-05 Thread Guillaume Drolet
Hi John,

I tried pprepair and imported the fixed shapefile into PostGIS. Although
pprepair said all polygons were fixed, there are still 22 slivers in my
topology. Maybe due to the value I used for tolerance? Finding the right
tolerance for building a topology isn't an easy task...



--
View this message in context: 
http://postgis.17.x6.nabble.com/Topology-cannot-delete-slivers-or-gaps-tp5007250p5007279.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


Re: [postgis-users] Topology: cannot delete slivers (or gaps)

2014-11-05 Thread Guillaume Drolet
Sandro,

Using your help, I filled most holes but three I couldn't fill using the
"UPDATE ... SET topogeom ..." approach. Will explore other approaches later.

Re the presence of holes, I may have deleted faces when running my function
for removing extra nodes. This needs checking.

Advancing in my work, I also got to drop edge 6341 and I end up with this:

https://dl.dropboxusercontent.com/u/5196336/edge6340.bmp

So the next step is to split edge 6340 about its middle, with the goal of
then connecting the newly created node with node 4753, thus creating an edge
where there should be one (i.e. on the light grey line). 

However, I get the error: *point not on edge* when running ST_ModEdgeSplit.
First, I tried to find the new node location using the coordinate capture
tool in QGIS but because I got that error I thought maybe I really need to
find a point that's on the edge so I used this approach:

/WITH edge AS (SELECT geom line
  FROM de_20k_topo.edge_data
  WHERE edge_id = 6340
), point AS (SELECT e.line, ST_ClosestPoint(e.line, 
ST_GeomFromText('POINT(-208737.346 811681.163)', 32198)) 
point
 FROM edge e
)
SELECT ST_ModEdgeSplit('de_20k_topo', 6340, point)
FROM point;

/

I still get the same error: point not on node.

How is one supposed to find a point that IS on an edge?



--
View this message in context: 
http://postgis.17.x6.nabble.com/Topology-cannot-delete-slivers-or-gaps-tp5007250p5007281.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


Re: [postgis-users] Topology: cannot delete slivers (or gaps)

2014-11-06 Thread Guillaume Drolet
Sandro: TopoGeo_addPoint sounds like the appropriate choice. I hadn't seen
this alternative topology function. Will try that, thanks.

Remi: Thanks for your advice on when to use what tools. I should give GRASS
a shot at topology cleaning.  Since topology-related work isn't something I
do everyday and PostGIS is a tool I use frequently, I naturally went in that
direction.



--
View this message in context: 
http://postgis.17.x6.nabble.com/Topology-cannot-delete-slivers-or-gaps-tp5007250p5007284.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


Re: [postgis-users] Topology: cannot delete slivers (or gaps)

2014-11-07 Thread Guillaume Drolet
I finally finished cleaning my topology, using a mix of SQL commands,
functions, and also using Sandro's QGIS topology editor plugin. 

In parallel, I also tried going the GRASS route, importing my original
shapefile in GRASS and applying the cleaning function (v.clean ...) as
suggested by Remi. I injected the result into PostGIS and created a
topology. Most of the 23 slivers were there in the topology so one of two
things: 1) GRASS didn't fix them (maybe I didn't choose an appropriate
tolerance) or 2) they were created by PostGIS when building the topology.

Anyway, I got what I finally got what I need, thanks a lot to all of you who
helped me.

I'm gonna ask you one last thing if I may: I want to replace the original
geom column (i.e. that used to build the topogeom) with a topologically
correct one, and keep all the associated ecological attributes. Is this the
right way to do it:

UPDATE syshiera.de_20k SET geom = topogeom::geometry; ?

Thanks.



--
View this message in context: 
http://postgis.17.x6.nabble.com/Topology-cannot-delete-slivers-or-gaps-tp5007250p5007286.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


[postgis-users] Dynamic parameters to ST_MapAlgebra

2014-12-10 Thread Guillaume Drolet
Hi,

I'm trying to find the best way to dynamically apply scaling factors to a
raster. I've tried different approaches using ST_MapAlgebra, both callback
and expression versions. The fastest method so far is this one but it
doesn't allow for dynamic parameters:

SELECT ST_MapAlgebra(rast, 1, '32BF', '([rast] * 0.01) + 0') AS rast,
FROM evi;

Compared to using the following function, the query above is ten times
faster:

CREATE OR REPLACE FUNCTION public.scale_sds_plv8(value double
precision[][][], 
pos integer[][], VARIADIC userargs text[])
RETURNS double precision AS
$$
  var g = Number(userargs[0]);
  var o = Number(userargs[1]);
  
  return (value * g + o);
   
$$ LANGUAGE plv8 IMMUTABLE;

WITH scaling_params AS (
SELECT gain, off_set
FROM meta
WHERE product = 'MCD12Q2' AND sds = 'EVI'
) 
SELECT ST_MapAlgebra(rast, 1, 'scale_sds_plv8(double
precision[], integer[], text[])'::regprocedure, '32BF', 'FIRST',
NULL::raster, 0, 0, 
VARIADIC ARRAY[gain, off_set]::text[]) AS rast 
FROM evi, scaling_params 

I want to take advantage of the speed of the first method above but be able
to pass parameters dynamically instead of hard-coding them in the query. To
this aim I tried this approach:

CREATE OR REPLACE FUNCTION scale_sds_plpgsql(rast raster, gain float8,
offs float8 DEFAULT 0)
RETURNS raster AS $$
DECLARE
sql text;
BEGIN
sql := 'SELECT ST_MapAlgebra(' || rast || 
', 1, ''32BF'', ''([rast] * ' ||
gain || ') + ' || offs || ''')';

EXECUTE sql;
END;
$$ LANGUAGE plpgsql IMMUTABLE;

SELECT scale_sds_plpgsql(rast, 0.01, 0.0)
FROM evi;

Running this last query, I get a synthax error which I don't understand:  

ERREUR: erreur de syntaxe sur ou près de «
F4000407FC05E10B2668C203F41A0BD86AA1F7C2
CONTEXT:  fonction PL/pgsql scale_rast(raster,double precision,double
precision), ligne 12 à instruction EXECUTE
** Error **

I'm sure there are experienced folks on this list who will have some clues
to solve this or who will come up with a better approach. 

Thanks a lot for your time.

G











--
View this message in context: 
http://postgis.17.x6.nabble.com/Dynamic-parameters-to-ST-MapAlgebra-tp5007468.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users

Re: [postgis-users] Dynamic parameters to ST_MapAlgebra

2014-12-10 Thread Guillaume Drolet
Hugues: thank a lot for your fix.

The function now runs without error... but it's been running for almost an
hour now :(

Compared to the 30-ish seconds it takes using the hard-coded query, and the
some 120 seconds using the PLV8 version, I cannot think of using the PLPGSQL
version in its current formulation. 

Does anyone else has a possibly faster suggestion, more in the ballpark of
the two other methods?

Thanks,

G



--
View this message in context: 
http://postgis.17.x6.nabble.com/Dynamic-parameters-to-ST-MapAlgebra-tp5007468p5007470.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users


Re: [postgis-users] Dynamic parameters to ST_MapAlgebra

2014-12-11 Thread Guillaume Drolet
Hi Regina,

"
Paragon Corporation-2 wrote
> Guillaume,
> 
> If the expression approach is faster, is there some reason you aren't just
> doing
> 
> 
> WITH scaling_params AS (
> SELECT gain::text As gain, off_set::text as off_set
> FROM meta
> WHERE product = 'MCD12Q2' AND sds = 'EVI'
> )
> SELECT ST_MapAlgebra(rast, 1, '32BF', '([rast] * ' || gain || ') + '
> ||
> off_set ) AS rast
> FROM evi, scaling_params;
> 
> Or is that too hard-coded for you?

No it isn't. My SQL skills are just not enough sharp yet. I can't say I have 
loads of experience with this language and some of the possibilities still
slip 
my mind. Thanks for the enlightening me!

> Now regarding the one that Hugues gave you, you have to run it differently
> from how you were calling your others, otherwise it would take forever
> since
> for every raster row you have it would apply map algebra to all the rows
> in
> your table.  So how did you call his?

Called it pretty much the same way as I showed in my first post: 

SELECT scale_sds_plpgsql(rast, 0.01, 0.0)
FROM evi
WHERE date = mydate;

> If you wanted to use the syntax you have I would rewrite your plpgsql as
> follows (and not use dynamic sql) - you don't need dynamic sql to do this
> just call ST_MapAlgebra directly.

When I referred to 'dynamic', I meant with passing parameters resulting from
a query. Maybe 'dynamic' isn't the right word?

Thanks for all your suggestions, I appreciate it. Will try them and see how
they improve things.

Best wishes,

G















--
View this message in context: 
http://postgis.17.x6.nabble.com/Dynamic-parameters-to-ST-MapAlgebra-tp5007468p5007475.html
Sent from the PostGIS - User mailing list archive at Nabble.com.
___
postgis-users mailing list
postgis-users@lists.osgeo.org
http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users