locate_along_measure: wrong values, invalid data ?

[l...@…]

I get absurd wrong results with locate_along_measure. Following up, force_2d crashes the server process.

Version Information: PostgresSQL 8.3.1 (binaries downloaded) on Windows XP Pro PostGIS 1.3.2, with Stack Builder

(I got the same behaviour with a older PostgresSQL 8.1)

---

---

create table t2 (id serial not null primary key); select AddGeometryColumn('public', 't2', 'geom', 31293, 'LINESTRINGM', 3);

insert into t2(geom) values (ST_GeomFromEWKT(' SRID=31293;LINESTRINGM( 6193.76 5337404.95 4519, 6220.13 5337367.145 4566, 6240.889 5337337.383 4603 )'));

select asewkt(st_locate_along_measure(geom, 4566)), * from t2;

---

---

select asewkt(force_2d(st_locate_along_measure(geom, 4566))), * from t2;

---

[pramsey]

Confirmed. Works for every measure EXCEPT 4566... right on the vertex. Seems to generally fail for values that are exactly interior vertices.

[pramsey]

This seems to be happening in lwcollection_construct, but it's not clear to me why. All the coordinates are correct up until that point. Ugly.

[mcayland]

Hmmmm. The problem appears to be that clip_seg_by_m_range() seems to think that the line segment between each 2 consecutive points requires clipping at the vertex point where the M values equals the measure we are searching for. Hence this returns the wrong information back to ptarray_locate_between_m() which then generates an invalid geometry.

Drilling down into clip_seg_by_m_range() I find that the problem is being caused by this piece of code here:

/*

• if m0 and m1 have the same value
• avoid computation of second point to
• reduce rouding problems with floating
• numbers. *
• The two points must be equal anyway. */

if ( m0 == m1 ) {

memcpy(p2, p1, sizeof(POINT4D)); if ( swapped ) ret |= 0x0010; else ret |= 0x0100; return ret;

}

Now AFAICT is totally bogus, since m0 and m1 is the range of measure being searched for - so why on earth should searching for an exact value make you want to make both points the same and update the clipping flags??!

Removing this piece of code completely seems to resolve the issue but I'm not sure why it was added in the first place. Paul, any ideas?

ATB,

Mark.

[pramsey]

When I traced this, the invalid values didn't show up until much farther along, so I am a little agog that this works, but it does.

Also, the logic of short circuiting doesn't seem bogus to me, but perhaps I am misreading... m0==m1 implies p2==p1, so why not?

If I leave this in place and simply print out p1 and p2 after this function runs, they are both valid (identical) values...

Actually, I'm obviously unclear on something, since the mechanism by which this returns points rather than lines isn't obvious to me.

[mcayland]

I think it becomes clearer once you actually look at the values that get passed into clip_seg_by_m_range(). Here are the values for the first segment from the example above:

p1->x=6193.76 p1->y=5.3374e+06 p1->m=4519

p2->x=6220.13 p2->y=5.33737e+06 p2->m=4566

m0=m1=4366 (the value we are searching for)

From the above, it's reasonably obvious that the M range we are searching to match (m0 to m1) should not have any destructive effect on the points it is testing at all, let alone override the clipping result.

I wonder if the original coder though the same as you at first glance, i.e. m0==m1 implies p2==p1? AFAICT the actual test for equality should be p1->m==p2->m implies p2==p1, but this case is already handled at the top of the function.

ATB,

Mark.

[mcayland]

Oops: of course the line above should read:

m0=m1=4566

ATB,

Mark.

[mcayland]

Hi Paul, since you are the owner of the bug, I'll leave the fate of this one up to you ;) Incidentally I located this fairly quickly using the new debugging infrastructure in SVN trunk, so you may like to try using trunk and running configure with --enable-debug to see the extra detail.

ATB,

Mark.

[l...@…]

the problem seems that with the m-value given and the way clip_seg_by_m_range() works two result elements (instead of one) are created.

i hope the attached fix is acceptable, it does the right thing for me.

Richard

[mcayland]

Thanks for the patch.

It looks like your solution and my solution are fairly similar, although my approach above involves removing the offending section of the code rather than skipping it with continue ;)

Does removing the section of code listed in comment 3 resolve the problem aswell? If so, I think that this is the preferred solution.

ATB,

Mark.

[l...@…]

Hi Mark,

my solution passes the regression tests, your solution modifies the results.

select 'LINEZM_6', asewkt(locate_between_measures('LINESTRING(0 10 10 40, 10 0 0 0)', 2, 2));

result as expected and with my solution LINEZM_6|POINT(9.5 0.5 0.5 2)

result with your solution LINEZM_6|LINESTRING(9.5 0.5 0.5 2,9.5 0.5 0.5 2)

i didn't find the specification of the semantics, so i can't decide what is the correct result. The POINT(...) result makes more sense for me.

However, I hope that any solution is making it into to the next release code, because a core dumping RDBMS isn't the thing i want rely on - i assume we can agree on this.

I attached a patch with some additional regression tests.

Kind regards, Richard

[mcayland]

Hi Richard,

I'm afraid I don't see that result here at all. Here are my results using SVN PostGIS 1.3 branch:

Current behaviour:

postgis13=# select asewkt(st_locate_along_measure(ST_GeomFromEWKT(' SRID=31293;LINESTRINGM( 6193.76 5337404.95 4519, 6220.13 5337367.145 4566, 6240.889 5337337.383 4603 )'), 4566));

asewkt

---

SRID=31293;MULTIPOINTM(1.99915902736236e+37 -1.05616313798637e+45

5.41519507362328e-315,) (1 row)

postgis13=# select 'LINEZM_6', asewkt(locate_between_measures('LINESTRING(0 10 10 40, 10 0 0 0)', 2, 2));

?column? | asewkt

---

LINEZM_6 | POINT(9.5 0.5 0.5 2)

(1 row)

Behaviour with attached patch:

postgis13=# select asewkt(st_locate_along_measure(ST_GeomFromEWKT(' SRID=31293;LINESTRINGM( 6193.76 5337404.95 4519, 6220.13 5337367.145 4566, 6240.889 5337337.383 4603 )'), 4566));

asewkt

---

SRID=31293;POINTM(6220.13 5337367.145 4566)

(1 row)

postgis13=# select 'LINEZM_6', asewkt(locate_between_measures('LINESTRING(0 10 10 40, 10 0 0 0)', 2, 2));

?column? | asewkt

---

LINEZM_6 | POINT(9.5 0.5 0.5 2)

(1 row)

Furthermore with the same patch I see the same results for your extra regression tests:

postgis13=# select 'LINEZM_7', asewkt(locate_between_measures('LINESTRING(0 0 0 0, 1 1 1 1, 2 2 2 2, 3 4 5 6)', 1, 2));

?column? | asewkt

---

LINEZM_7 | LINESTRING(1 1 1 1,2 2 2 2)

(1 row)

postgis13=# select 'LINEZM_8', asewkt(locate_between_measures('LINESTRING(0 0 0 0, 1 1 1 1, 2 2 2 2, 3 4 5 6)', 1.1, 2.1));

?column? | asewkt

---

LINEZM_8 | LINESTRING(1.1 1.1 1.1 1.1,2 2 2 2,2.025 2.05 2.075 2.1)

(1 row)

postgis13=# select 'LINEZM_9', asewkt(locate_between_measures('LINESTRING(0 0 0 0, 1 1 1 1, 2 2 2 2, 3 4 5 6)', 2, 2));

?column? | asewkt

---

LINEZM_9 | POINT(2 2 2 2)

(1 row)

postgis13=# select 'LINEPZM_1', asewkt(locate_along_measure('LINESTRING(0 0 0 0, 1 1 1 1, 2 2 2 2, 3 4 5 6)', 0));

?column? | asewkt

---

LINEPZM_1 | POINT(0 0 0 0)

(1 row)

postgis13=# select 'LINEPZM_2', asewkt(locate_along_measure('LINESTRING(0 0 0 0, 1 1 1 1, 2 2 2 2, 3 4 5 6)', 1));

?column? | asewkt

---

LINEPZM_2 | POINT(1 1 1 1)

(1 row)

postgis13=# select 'LINEPZM_3', asewkt(locate_along_measure('LINESTRING(0 0 0 0, 1 1 1 1, 2 2 2 2, 3 4 5 6)', 1.5));

?column? | asewkt

---

LINEPZM_3 | POINT(1.5 1.5 1.5 1.5)

(1 row)

So I'm really confused as to why your results don't match mine?!

ATB,

Mark.

[l...@…]

Hi Mark,

you are right. I assume that I messed up my codebase a little.

I can reproduce your results with a clean svn checkout and with the postgis-1.3.3.tar.gz release package.

Kind regards, Richard

[mcayland]

Thanks for testing again Richard. Are you happy for me to commit this fix to 1.3 and trunk?

ATB,

Mark.

[mcayland]

Hi Richard,

I've applied my patch from above to both the 1.3 and trunk versions, and included your enhanced regression tests to monitor these conditions in future. Thanks for your help with this bug report.

Many thanks,

Mark.

[l...@…]

​http://postgis.googlecode.com/issues/attachment?aid=-759732332149041948&name=issue-21-patch.diff

[mcayland]

I know that this is a old bug (and even marked as closed), but I've been investigating Stephen Davies' regression failure and have found that on some platforms, removing the logic in comment #3 re-introduces floating point errors :(

As a result of this, I've realised that the combination of memcpy/clipping logic was wrong and so have spent quite a bit of time testing this on Stephen's machine. With the latest SVN commits to 1.3 branch and trunk, I've re-introduced the memcpy() and adjusted the logic so that his machine now passes regression.

However: the LRS logic is complicated. If you have some LRS test cases lying around, please run them through a recent SVN update and make sure I haven't missed any corner cases.

I don't see this as an emergency fix to push out a 1.3.5 yet, but I've backpatched it as it is likely to catch someone out in the future.

ATB,

Mark.

[l...@…]

​http://postgis.googlecode.com/issues/attachment?aid=-140617382247917647&name=regression_issue-21-patch.diff

[l...@…]

​http://postgis.googlecode.com/issues/attachment?aid=4989425466285116606&name=lwgeom_functions_lrs.patch