[postgis-tickets] r17137 - Change point-in-cone routine to be more tolerant
Paul Ramsey
pramsey at cleverelephant.ca
Fri Jan 11 01:07:17 PST 2019
Author: pramsey
Date: 2019-01-11 13:07:17 -0800 (Fri, 11 Jan 2019)
New Revision: 17137
Modified:
trunk/liblwgeom/cunit/cu_tree.c
trunk/liblwgeom/lwgeodetic.c
trunk/liblwgeom/lwgeodetic_tree.c
trunk/liblwgeom/lwgeodetic_tree.h
Log:
Change point-in-cone routine to be more tolerant
of narrow edges, and hopefully more robust.
References #4290
Modified: trunk/liblwgeom/cunit/cu_tree.c
===================================================================
--- trunk/liblwgeom/cunit/cu_tree.c 2019-01-11 13:50:49 UTC (rev 17136)
+++ trunk/liblwgeom/cunit/cu_tree.c 2019-01-11 21:07:17 UTC (rev 17137)
@@ -59,13 +59,13 @@
/* Point in square */
g = lwgeom_as_lwline(lwgeom_from_wkt("LINESTRING(-1 -1,1 -1,1 1,-1 1,-1 -1)", LW_PARSER_CHECK_NONE));
c = circ_tree_new(g->points);
- rv = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
CU_ASSERT_EQUAL(rv, 1);
/* Point on other side of square */
pt.x = 2.0;
pt.y = 0.0;
- rv = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
CU_ASSERT_EQUAL(rv, 0);
/* Clean and do new shape */
@@ -78,13 +78,13 @@
g = lwgeom_as_lwline(lwgeom_from_wkt("LINESTRING(-1 -1,0 -1,1 -1,1 0,1 1,0 1,-1 1,-1 0,-1 -1)", LW_PARSER_CHECK_NONE));
c = circ_tree_new(g->points);
//circ_tree_print(c, 0);
- rv = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
CU_ASSERT_EQUAL(rv, 1);
/* Point on other side of square, stab passing through vertex */
pt.x = 2.0;
pt.y = 0.0;
- rv = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
CU_ASSERT_EQUAL(rv, 0);
/* Clean and do new shape */
@@ -97,7 +97,7 @@
g = lwgeom_as_lwline(lwgeom_from_wkt("LINESTRING(-1 -1,0 -1,1 -1,1 0,1 1,0 0,-1 1,-1 0,-1 -1)", LW_PARSER_CHECK_NONE));
c = circ_tree_new(g->points);
//circ_tree_print(c, 0);
- rv = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
//printf("rv %d\n", rv);
CU_ASSERT_EQUAL(rv, 0);
@@ -104,7 +104,7 @@
/* Point inside "w" thing, stab passing through vertexes and grazing pointy thing */
pt.x = 0.8;
pt.y = 0.0;
- rv = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
//printf("rv %d\n", rv);
CU_ASSERT_EQUAL(rv, 1);
@@ -134,7 +134,7 @@
c = circ_tree_new(p->rings[0]);
//circ_tree_print(c, 0);
rv_classic = lwpoly_covers_point2d(p, &pt);
- rv_tree = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv_tree = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
CU_ASSERT_EQUAL(rv_tree, rv_classic);
circ_tree_free(c);
lwgeom_free(g);
@@ -148,7 +148,7 @@
//printf("OUTSIDE POINT(%g %g)\n", pt_outside.x, pt_outside.y);
c = circ_tree_new(p->rings[0]);
rv_classic = lwpoly_covers_point2d(p, &pt);
- rv_tree = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv_tree = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
CU_ASSERT_EQUAL(rv_tree, rv_classic);
circ_tree_free(c);
lwpoint_free(lwpt);
@@ -163,7 +163,7 @@
//printf("OUTSIDE POINT(%g %g)\n", pt_outside.x, pt_outside.y);
c = circ_tree_new(p->rings[0]);
rv_classic = lwpoly_covers_point2d(p, &pt);
- rv_tree = circ_tree_contains_point(c, &pt, &pt_outside, &on_boundary);
+ rv_tree = circ_tree_contains_point(c, &pt, &pt_outside, 0, &on_boundary);
CU_ASSERT_EQUAL(rv_tree, rv_classic);
circ_tree_free(c);
lwpoint_free(lwpt);
@@ -349,6 +349,23 @@
CU_ASSERT_DOUBLE_EQUAL(d1, d3, 0.00000001);
CU_ASSERT_DOUBLE_EQUAL(d1, d4, 0.00000001);
+
+ /* Ticket #4290 */
+ const char *wkb1 = "0103000020E610000001000000EB010000560E7951F83853C0986D8DC811AE42404601F3C2F93853C028F9A6C80FAE42400EF46097FB3853C07B5CA9550FAE4240CC9BF16EFD3853C0F0D49BDC0FAE42400662BD42FF3853C03DF79F370FAE4240FC814915013953C0669AAA2E0EAE42401B6D54A7023953C0372CBC910CAE4240A451EF59043953C0336BCB570BAE42409C75E90B063953C08A26DAEB09AE42401554B6DF073953C04131DD4609AE4240F900D2D2093953C006A1DFA008AE424074706C850B3953C09198EB6607AE42401D7798580D3953C06B9AF28F06AE4240D75FC42B0F3953C0B5E9F6B805AE42409732F0FE103953C06600FBE104AE4240A0475DD3123953C0198EFA6E04AE4240B9805AC8143953C0AA9EF25E04AE42407E92B8BC163953C0E888EC1C04AE4240DA1615B1183953C0FBB7E8DA03AE4240184CE1841A3953C0B87CE93503AE42402370AD581C3953C07508EA9002AE4240B1483C301E3953C0247BD81703AE42405ECEE5321F3953C02740AB9905AE424008B91F18203953C0691E78B208AE42402245F374213953C07DE45A050AAE42408017844C233953C0E4A0488C0AAE4240005FE240253953C08116434A0AAE424029093F35273953C011DD3A080AAE4240DA0AAC09293953C0FBAA379509AE42403CA208FE2A3953
C041F22E5309AE424033BED4D12C3953C090742DAE08AE4240E30392C52E3953C080D5283A08AE424020960D7A303953C080A2279607AE424068FB69AC323953C0F3911E2007AE4240A3059541343953C0F3B51D7D06AE4240B2CF6015363953C0FD161BD805AE424027B38BE8373953C079E51C0105AE424051CB55BC393953C06CE1195C04AE4240397D808F3B3953C08DC4188503AE424014BD5A433D3953C09BC41AAF02AE4240F98D76363F3953C0A52D160902AE4240A0E9410A413953C04739126401AE4240BA535A3B433953C0339C0D8A00AE4240000C95B0443953C0F58710B6FFAD4240C73EBF83463953C010170EDFFEAD4240EDDCF774483953C03D5111A3FDAD42400C5120484A3953C06C780ECCFCAD4240B8F7B7FA4B3953C0012C1392FBAD424011D740CD4D3953C046821489FAAD4240AF9BC99F4F3953C055A01580F9AD4240443CF372513953C0FBDC11A9F8AD4240D6CE1C46533953C0F95A10D2F7AD42405CB64419553953C026330CFBF6AD42402417CDEB563953C0C56F0CF2F5AD42400055B4BD583953C032A40EB7F4AD4240191E4B705A3953C078B9117DF3AD42401FF250025C3953C024881CE0F1AD4240131655945D3953C0813A2743F0AD424094704A465F3953C0C4F62BD7EEAD4240321FD218613953C0EDFA2ACEEDAD424017AB59EB623953C0EA4
C27C5ECAD42403B24E1BD643953C0A6E025BCEBAD424081351871663953C0DC6D25B4EAAD4240EA11AE23683953C0DBFB267AE9AD4240743A35F6693953C0C5EE2471E8AD424020615DC96B3953C083FD1D9AE7AD424054EAF27B6D3953C06CF21E60E6AD42408606BC4F6F3953C0D36115BBE5AD4240183BA2E3703953C0AE0615B4E4AD4240B2EA1A98723953C048570C10E4AD42402E54006A743953C0583A0BD5E2AD4240C559951C763953C008380B9BE1AD4240F49F290D783953C038060D2DE0AD424025BF2B9F793953C0BD5C1290DEAD424002A610717B3953C0A7681055DDAD4240B1FE03237D3953C09E4B14E9DBAD42405E6798D57E3953C0805013AFDAAD42407DC73C68803953C033B41544D9AD42409343ED18823953C0BE721D74D7AD4240E7634DAA833953C066EA28A5D5AD4240DF7CCDFB843953C000F23874D3AD424085E94B4D863953C0F9634B43D1AD4240D6D3E85E873953C0017D67B0CEAD424065B12733883953C008167D83CCAD42406A881464893953C0B39497EFC9AD42409C16F1B48A3953C0B2DDAB8CC7AD42404B7ACD058C3953C0930AC029C5AD42402AFF4A578D3953C011E6D1F8C2AD42402D3037888E3953C0E201EC64C0AD42403D5DB4D98F3953C0B12DFB33BEAD42406D1950EB903953C0C8A216A1BBAD42400D889BDD913953C00C17310F
B9AD4240E773B570923953C0D675591CB6AD4240F9CED003933953C0C0487F29B3AD4240B6B2CA14943953C028BE9C64B0AD42406EB34928953953C04056AC67AEAD42405F902479963953C0C7A5BF04ACAD4240522FAFAA973953C02B70D4A2A9AD4240BA954C7E983953C0A04EF043A7AD42404764A48E993953C0C59C0F4DA4AD4240B2ABEE809A3953C0228F29BBA1AD4240A10049539B3953C0FEB449F89EAD424080AAA1259C3953C0FEDC69359CAD4240321C4DD99C3953C0755889A599AD424030BB05AB9D3953C07396ABB096AD42409B990F5E9E3953C0A237CDEE93AD4240D40FD50F9F3953C0E64AF3C890AD4240324B90A39F3953C0B17616088EAD424032895817A03953C000883D168BAD4240B4FFB4ABA03953C0FD735E8788AD4240D4CD0F40A13953C08FED7CF885AD424060BB77F2A13953C0439EA00483AD42400D96DFA4A23953C053C1C61080AD424023602897A33953C0BC2CE07E7DAD4240A1AACF88A43953C09E3AFEBA7AAD4240068DD83BA53953C0C1901FF977AD42402556E1EEA53953C03DDF403775AD4240AB60AA62A63953C0A8B7674572AD42405D96D397A63953C0F0F491556FAD42400A935E8EA63953C0238CBF676CAD4240C9B6B671A63953C029C16D3E6BAD4240C2815AF3A53953C04C461B6269AD4240F7E145ABA53953C02EB84E7666AD4
240878F2FA1A53953C0FE817E5663AD4240D6525AD6A53953C0F3B1A86660AD424046A3214AA63953C03E10D2745DAD4240E367EABDA63953C095E4F8825AAD42408BD54090A73953C0957B18C057AD4240353D5923A83953C036913DCD54AD4240AC72AFF5A83953C045165D0A52AD4240A942F8E7A93953C0261A76784FAD424012D10E7BAA3953C09C279B854CAD4240FBAE644DAB3953C06F8EBAC249AD424037CDCAFFAB3953C027D2DDCE46AD42405D2D4354AC3953C07AB608DE43AD42409529C7FFAC3953C075B49F9D41AD4240C07BAD14AD3953C0F90522123FAD42405A18380BAD3953C0369C4F243CAD424045502101AD3953C035E2810439AD42401DA0E774AD3953C08793A81236AD42408EAEEE27AE3953C05D89C95033AD4240401206BBAE3953C0C263EE5D30AD4240E30A7E0FAF3953C0EEBE166D2DAD4240C2C293A2AF3953C0E2173E7A2AAD4240C45CA935B03953C0BCF1628727AD4240B0F7AFE8B03953C09EC983C524AD42401864C57BB13953C01098A8D221AD4240D1AF99CCB23953C05FA5B96F1FAD4240026AF15DB43953C0C7DDBDA01DAD424016EAE8EFB53953C017C3BF031CAD4240A85E74A2B73953C022B6B8C91AAD42407E110D35B93953C06694B55E19AD42401F4264C6BA3953C0C19FBB8F17AD4240C39EC837BC3953C0A42AC38F15AD424003
B05C49BD3953C0753BD9FC12AD42407E24443CBE3953C0EB50F19C10AD42404E4949EFBE3953C055140FDB0DAD4240AC5C9CC1BF3953C078722D180BAD424096DC3EB3C03953C0C2DF495408AD4240BA3E9385C13953C0491A689105AD4240FD8A2697C23953C00F4B80FE02AD4240070823C6C33953C0DB329DD4FFAC42407B0EC979C43953C02BFBBA44FDAC42400D54C96AC53953C05F5AD94EFAAC42407A97CD1DC63953C00249F98CF7AC42409C3EF190C63953C06D782169F4AC42406D7417C6C63953C00E514B79F1AC42407394DC39C73953C021BD7387EEAC4240A380EFCCC73953C0D4C49794EBAC42408BC6429FC83953C010A4B5D1E8AC4240ED2FB5F0C93953C04C15C0A0E6AC4240AC823622CB3953C0BE91D33EE4AC42401756FC54CC3953C0F384DD40E2AC424072E300C7CD3953C0B034E172E0AC42403982D95ACF3953C0176DD76BDFAC4240656AF22DD13953C022BAC494DEAC424097FB45BFD23953C0CAEAC5C5DCAC4240CEA409F2D33953C085CED1C7DAAC424024871C44D53953C04F11D9C8D8AC424052FB5A74D63953C00718F002D6AC4240E0CF4C47D73953C0B20A0B72D3AC424095902C78D83953C065371DDED0AC4240BFDB200ADA3953C00A731B41CFAC4240648B95DCDB3953C0B73C0C38CEAC42406550ADAFDD3953C0060DF860CDAC42408C03D46
2DF3953C0107EEA58CCAC42406D1C35D4E03953C01BA1EE58CAAC42407F610325E23953C007BCFBF5C7AC424070453AB8E33953C06584F4BCC6AC4240C7C9BE6AE53953C00685E882C5AC4240712F343DE73953C01F6BD579C4AC4240143C4B10E93953C05F57C2A2C3AC4240C35902E4EA3953C01462AAFDC2AC424049F3BAB7EC3953C004279258C2AC42407834158CEE3953C0D4FC74E5C1AC4240F8B5CD5FF03953C090C95E40C1AC42409158E432F23953C000204869C0AC42401B1368E5F33953C0FB803A2FBFAC42407AB0EB97F53953C081B12CF5BDAC4240B3D5A36BF73953C07C291350BDAC4240D68DB83EF93953C0ACB2FB78BCAC424097124734FB3953C0FBE9D39ABCAC424018E6A6F2FB3953C08955BA5CBDAC424051AB0417FC3953C054999EEBBEAC4240ACDE61A4FA3953C00A82A187C0AC42406D51DEF1F83953C005E5AFC1C1AC424031F8EB5FF73953C009D9B45EC3AC42400C11FBCDF53953C03111BCFBC4AC4240DF71083CF43953C033B3C098C6AC4240AC24268AF23953C045A9C904C8AC424067D6E016F13953C09862D06EC9AC424099E43286EF3953C07703D06FCBAC4240FF1B92D5ED3953C07C8ED13FCDAC42403D7F6104EC3953C0477CDEACCEAC4240EDD2B8CFEA3953C0752BE014D0AC424026866B00E93953C0B322E117D2AC42403B26C78DE739
53C02740E2B3D3AC42404D3B241BE63953C0FEA7E54FD5AC4240F4E8D189E43953C0ACFAE31ED7AC42401C5ADFF7E23953C00147E9BBD8AC4240E400FA45E13953C09E61F027DAAC4240301B1694DF3953C047B9F993DBAC4240851C2302DE3953C0D109FC30DDAC4240317E1F90DC3953C02CE4F9FEDEAC424032055F1FDB3953C021ADF030E1AC42400B3F0D8ED93953C0694FF0FFE2AC4240E8D4AA1CD83953C06185EBFFE4AC4240583B48ABD63953C0301FE4FFE6AC42407272F419D53953C0E85AE3CEE8AC4240BC0AD2E7D33953C08758D5FEEAAC42408D904FB6D23953C00112C560EDAC4240117DCE84D13953C0B9A7B4C2EFAC424032B64B53D03953C03733A424F2AC4240A343D901CF3953C03AD89755F4AC424089BC57D0CD3953C0712687B7F6AC4240CA6EE37ECC3953C066A27AE8F8AC42403019516DCB3953C0791D637BFBAC42408A147D5ACA3953C074DF54AAFDAC42408A942D4AC93953C0054A36A100AD4240C7F7B9E7C73953C08FDE2D6004AD42403B4F27CEC63953C0274EF85907AD4240AF45811AC63953C075A0DAE909AD4240EB958F09C53953C0FD75C0AE0CAD424048D14755C43953C05D67A70C0FAD42408667F744C33953C063E68A0312AD4240A3DE5453C23953C0B1A16EC714AD424018496080C13953C0A018555817AD424070E1788DC03953C00
B373DB819AD42402DCE965CBF3953C09D38294C1CAD4240FA99024BBE3953C007D410DF1EAD4240A4257D19BD3953C03EEFFE4021AD4240389A9AE8BB3953C0B2AEEAD423AD4240B6D714B7BA3953C0059CD83626AD42401C0C80A5B93953C0CEE7BFC928AD42402A42DCB3B83953C0B012A38D2BAD4240258D9A83B73953C088CC89532EAD42403F25A3B0B63953C0E3C26FE430AD4240CA891ABDB53953C019E25B1233AD42403E3776CBB43953C079D23ED635AD4240009C6F18B43953C045261E9838AD42407DE76865B33953C06F72FD593BAD4240B42162B2B23953C0E130DF1B3EAD4240BE634C1FB23953C081F2B70E41AD4240D60A358CB13953C0B135930144AD4240455E1019B13953C08E676A2547AD4240CCB305A4B03953C0EF4848B349AD4240F6ABEF10B03953C0057223A64CAD424025C8279DAF3953C067DDFC974FAD4240D89CAF48AF3953C0CB80D48852AD4240985D37F4AE3953C05DA8A97955AD4240BBE80EBFAE3953C023847F6958AD4240099C966AAE3953C07B22575A5BAD42401D70CEF6AD3953C0BB81304C5EAD42405BC00783AD3953C0F4D0093E61AD42405C4F0BD5AC3953C026ADCF8F65AD424083E84261AC3953C00101A98168AD4240D5FD7BEDAB3953C0DA4482736BAD42400DE66198AB3953C0560B5C326EAD42403516FBE5AA3953C0FCC338
2671AD42402B09E452AA3953C0553F111974AD424052491BDFA93953C04881EA0A77AD42405FA50F6AA93953C083A8CA9879AD4240DBB6E9F6A83953C0F2259FBC7CAD424059D5D063A83953C098167AAF7FAD4240F55C09F0A73953C0F13D53A182AD42409746407CA73953C0CB6E2C9385AD424087AD28E9A63953C06ECB048688AD4240E0725F75A63953C0F4F4DD778BAD4240B44EA6E1A53953C089F9BA388EAD4240DD863B6DA53953C0F45296F890AD4240794332BAA43953C0D47977BA93AD4240E26B2A07A43953C05712567C96AD424003CD1074A33953C0ACD6306F99AD4240563707C1A23953C00A6E0F319CAD4240849CBDCEA13953C03E80F3C29EAD42406A5165FCA03953C027A3D385A1AD4240FEFF6929A03953C01779B816A4AD42404869C1379F3953C0902D9ADAA6AD4240E0AF18469E3953C0F0D37B9EA9AD4240D2EECC539D3953C0A9296230ACAD42406E9282619C3953C071EA45C2AEAD4240991D386F9B3953C0F1162C54B1AD424086649D5D9A3953C06C5811E7B3AD4240A3AC526B993953C09366F778B6AD4240EE904B3C983953C07501D5A2B9AD424039DFBD48973953C0325CBFD0BBAD4240AA002075963953C0D7E9A52FBEAD4240DAC4D482953953C02BBD8BC1C0AD42409E608990943953C053086F53C3AD42408F39FC5E933953C0D3135AB5C5A
D4240F0F0614D923953C0EFBE3E48C8AD42402AAD243B913953C0561428A9CAAD4240B1159709903953C0BDDF120BCDAD42409D89F41F8E3953C02157F49ED0AD4240DD1885AE8C3953C04248E69ED2AD4240A90A075D8B3953C021DDD6CFD4AD42402740870B8A3953C079E8C400D7AD42401E22AABA883953C09E95B063D9AD42407056DB49873953C039D29F95DBAD4240DA6ECBD7853953C0C2C39563DDAD42408B0F5B66843953C0DEF08663DFAD4240530EF9D4823953C0AD177E32E1AD4240E62AE762813953C095AB7300E3AD4240012436B27F3953C054DE6BD0E4AD42408ABC823F7E3953C0625B636CE6AD42403E73E16E7C3953C066BD5E0BE8AD42401565EBFA7A3953C06CCD5C43E9AD4240E1BD6428793953C0AC67634CEAAD4240F5F3DD55773953C0C44F6755EBAD424055E1B582753953C0DA326D2CECAD42403B22EEAE733953C08F7D77D1ECAD424001EEC5DB713953C033697FA8EDAD42409CC23E09703953C0AE6F82B1EEAD42409C3B75356E3953C0B91C8C56EFAD42406CC54C626C3953C08F5E932DF0AD4240C48B76706A3953C08CE99637F1AD424006E24BDB683953C099879DDAF1AD4240E1CAC508673953C01D6B9FE3F2AD4240E2133E36653953C01D9DA3ECF3AD424048E46544633953C06748A6F6F4AD4240959C7EB0613953C0B33CA5FDF5AD4240
A28455DD5F3953C0C0FAAAD4F6AD42405C478B095E3953C0E3B0B279F7AD4240D9FD61365C3953C09A83B550F8AD4240FE3DDB635A3953C0A15AB859F9AD4240CE9C42B1583953C0DFD9B293FAAD42408473ABFE563953C0D595AFCDFBAD42409AD3222C553953C0F9D8B1D6FCAD424009718B79533953C03C30AE10FEAD4240631F04A7513953C04780AD19FFAD424070F66AF44F3953C0927FA95300AE4240F16DE3214E3953C02E63A85C01AE42409D3D5A4F4C3953C0511BA76502AE424086FAD07C4A3953C02D15A86E03AE424029A1A6A9483953C04C8BA84504AE42403014BED7463953C05D6CA48005AE4240A1B844E5443953C0CCCAA75806AE42405759B950433953C0F8CFA52D07AE4240A2FE7F9D413953C062CFA43508AE4240049BA5E93F3953C07E51A30B09AE4240D9D2CCB93D3953C03CB0A0490AAE424002DCD4453C3953C0DFB498810BAE4240028599923A3953C0BD7C94890CAE4240638100E0383953C045DA8DC30DAE4240B98BC62C373953C0D7AD8BCB0EAE42405623ED3A353953C07D6088D50FAE4240078D1087333953C0358A85AB10AE4240747735D3313953C0F0EF848111AE4240FD4F79DF2F3953C078EF89F511AE42401CD1AC0B2E3953C0F5D28B9A12AE42400741E0372C3953C0707D8D3F13AE42400A5106462A3953C0FD3F8B4914AE4240B75AD
A72283953C019D0872015AE42402B8F0D9F263953C051CC88C515AE42403977E1CB243953C0F464879C16AE424070CEB6F8223953C0523E837317AE4240F1ED1A46213953C012037AAD18AE42406B3B21941F3953C094606E191AAE4240A6EFE5E01D3953C0A6E166211BAE4240F9C8FA0E1C3953C0A09E5D5C1CAE4240B37CB07B1A3953C03FFE52951DAE42406D8864E8183953C0E14048CE1EAE42402D401B17173953C0BCB5393B20AE424051917E64153953C027272F7521AE4240130184B2133953C0923222E122AE4240CAA2E8FF113953C07636171B24AE424066F8DE6D103953C004C104B825AE4240B777E2BB0E3953C0324CF72327AE4240EA87D8290D3953C03E81E4C028AE424044AC3C770B3953C018CBD8FA29AE4240F6EB30E5093953C02EB7C5972BAE42408F0AD633083953C0ABDAB2352DAE424002FA2AA1063953C02606A4A02EAE42404ED11E0F053953C0AE73903D30AE4240BB12147D033953C074AB7CDA31AE4240EB3BA8EB013953C0DB9F66A933AE424039CB3D5A003953C0C95E507835AE42406EA6D1C8FE3853C0AD013A4737AE4240E355C636FD3853C0B69525E438AE424068D4C8C2FB3853C01581171C3AAE42403E445C31FA3853C05CAC00EB3BAE4240311AF19FF83853C035A2E9B93DAE424010B2E30DF73853C0FAA1D4563FAE42400FB4D77BF5
3853C0006CBFF340AE4240ED8B6AEAF33853C0FEF3A7C242AE4240DF7DDF98F23853C0100E89F344AE42406A161408F13853C03A196FF446AE42400894D897EF3853C0150A4E5849AE4240C6140C45EE3853C0508B35254BAE424044DA7EF3EC3853C0573816564DAE42408A6D0382EB3853C07F06F9554FAE4240652F95F0E93853C07E8CE02451AE4240BF9B689EE83853C0F26DC52353AE4240D910FA0CE73853C012A9ACF254AE424046087E9BE53853C0F4EA8EF256AE424007586FCBE33853C08F0F76C358AE4240BD93E279E23853C06FCD55F45AAE4240FED0A447E13853C0065F35245DAE424055620854E03853C0291E17525FAE424007307B02DF3853C0BC91F68261AE424051D17DD1DD3853C0BA1ECF1664AE424083A630BFDC3853C0FD1EAC7766AE42402981D4CCDB3853C0B2D7860969AE42400C7127BBDA3853C0DD175F9C6BAE42409A08CBC8D93853C030B2392E6EAE4240F7B01DB7D83853C0EECF11C170AE42406A3C70A5D73853C05755EC5373AE42400C9BC293D63853C0834EC4E675AE4240883C2562D53853C0A7CAA04878AE4240E892C66FD43853C0EDA678DA7AAE42404EEA283ED33853C0F0F9543C7DAE4240368E890CD23853C09C42319E7FAE424030F38BDBD03853C0D9440B3282AE4240DF71DDC9CF3853C09AC9E2C484AE424098FB2E7ACE3853C
00EFBB78B87AE4240808E3F67CD3853C0CB8F98BA89AE424017254174CC3853C00009721A8CAE4240436C3263CB3853C0EF9949DF8EAE4240057FBDCCC93853C0ACED431E8FAE424091FA48EBC83853C0207A62318DAE42406547459BC93853C00C6D9A7589AE42403DFF116DCA3853C01975C78086AE4240C714705FCB3853C0733EF0EE83AE424008E42D13CC3853C06956165F81AE42406690DC24CD3853C083133FCC7EAE424023A4BC95CE3853C0DA29619A7CAE42405ABD4BE7CF3853C0961783697AAE4240D8F7BC78D13853C051319E9A78AE42401D3FDDEAD23853C0EF21B9CC76AE424023BB7A1CD43853C0F6D6DC6A74AE4240981CE8EED43853C0392707A871AE424093F044E1D53853C054D72C166FAE4240EACEF2B4D63853C0AA5050B76CAE4240D0495E87D73853C0EF8B7AF469AE424026A6A8DAD73853C072ADAD9F66AE4240C0F551B1D73853C056C2DE8063AE4240C3526E6CD73853C0C891FB8E61AE424023C21A02D83853C0E9E71C645FAE4240D2669AF7D83853C0693535CC5DAE4240EE32DAE8DA3853C0D77446905CAE42405C718A7BDC3853C04FB759255BAE424010B966ECDD3853C0C2747AF358AE4240899E13C0DE3853C085849D9456AE42405B9A7F92DF3853C0A24AC7D153AE4240BD371805E03853C0C338FA7B50AE4240E8A511FBDF3853C0DC4D
2B5C4DAE42408987EC6EE03853C0E66A576A4AAE424015975641E13853C0732881A747AE4240071EB333E23853C09513A61545AE42402167AE26E33853C06AD9C8B542AE42400D9A28D9E33853C0C8ABF4C13FAE42406423438CE43853C0C1521E003DAE42406E935D3FE53853C00FF2473E3AAE4240125D76F2E53853C06396717C37AE424015C98003E73853C008BF9AB734AE42406E519E37E83853C0FBDFB11D33AE424055715AE8E93853C08D74CA4D31AE4240E6CC853AEB3853C0BA81E54E2FAE4240C97F518DEC3853C0DB4DFE812DAE42405BADDDDEED3853C096B31D512BAE4240C39418D3EE3853C0B07E3B5529AE42400AD53307F03853C0D2974FBB27AE4240C3950F3AF13853C00ACA6CBD25AE42408AA8B9CFF13853C022328D9223AE4240B1B39243F23853C0F7D3B8A020AE42400F123C55F33853C0BC8ADC0D1EAE424095D59447F43853C051EB027C1BAE4240F05D0FB9F53853C036D51C7C19AE4240FC4948EBF63853C075CB3B4C17AE4240CEB8609EF73853C05AA0648A14AE4240560E7951F83853C0986D8DC811AE4240";
+ const char *wkb2 = "0103000020E610000001000000080000002CAFE119EF2F53C02D6D4097DFBE42400050FBE7F32F53C05A08B32BD9BE42402ED0C3E9E52F53C0132FFFAE9DBE424048885052FE2F53C0DE5D06256BBE42409534CFEF163053C07D170F98ADBE42402E72095F303053C025ECD50CEABE424059CE3A7C103053C0235393C114BF42402CAFE119EF2F53C02D6D4097DFBE4240";
+ spheroid_init(&s, WGS84_MAJOR_AXIS, WGS84_MINOR_AXIS);
+ lwg1 = lwgeom_from_hexwkb(wkb1, LW_PARSER_CHECK_NONE);
+ lwg2 = lwgeom_from_hexwkb(wkb2, LW_PARSER_CHECK_NONE);
+ c1 = lwgeom_calculate_circ_tree(lwg1);
+ c2 = lwgeom_calculate_circ_tree(lwg2);
+ d1 = lwgeom_distance_spheroid(lwg1, lwg2, &s, threshold);
+ d2 = circ_tree_distance_tree(c1, c2, &s, threshold);
+ circ_tree_free(c1);
+ circ_tree_free(c2);
+ lwgeom_free(lwg1);
+ lwgeom_free(lwg2);
+ CU_ASSERT_DOUBLE_EQUAL(d1, 18358.866, 0.01);
+ CU_ASSERT_DOUBLE_EQUAL(d2, 18358.866, 0.01);
}
static void test_tree_circ_distance_threshold(void)
Modified: trunk/liblwgeom/lwgeodetic.c
===================================================================
--- trunk/liblwgeom/lwgeodetic.c 2019-01-11 13:50:49 UTC (rev 17136)
+++ trunk/liblwgeom/lwgeodetic.c 2019-01-11 21:07:17 UTC (rev 17137)
@@ -3380,13 +3380,46 @@
/* The projection of start onto the center defines the minimum similarity */
min_similarity = dot_product(A1, &AC);
- /* The projection of candidate p onto the center */
- similarity = dot_product(P, &AC);
+ /* If the edge is sufficiently curved, use the dot product test */
+ if (fabs(1.0 - min_similarity) > 1e-10 )
+ {
+ /* The projection of candidate p onto the center */
+ similarity = dot_product(P, &AC);
- /* If the point is more similar than the end, the point is in the cone */
- if ( similarity > min_similarity || fabs(similarity - min_similarity) < 2e-16 )
+ /* If the projection of the candidate is larger than */
+ /* the projection of the start point, the candidate */
+ /* must be closer to the center than the start, so */
+ /* therefor inside the cone */
+ if ( similarity > min_similarity )
+ {
+ return LW_TRUE;
+ }
+ else
+ {
+ return LW_FALSE;
+ }
+ }
+ else
{
- return LW_TRUE;
+ /* Where the edge is very narrow, the dot product test */
+ /* fails, but we can use the almost-planar nature of the */
+ /* problem space then to test if the vector from the */
+ /* candidate to the start point in a different direction */
+ /* to the vector from candidate to end point */
+ /* If so, then candidate is between start and end */
+ POINT3D PA1, PA2;
+ vector_difference(P, A1, &PA1);
+ vector_difference(P, A2, &PA2);
+ normalize(&PA1);
+ normalize(&PA2);
+ if (dot_product(&PA1, &PA2) < 0.0)
+ {
+ return LW_TRUE;
+ }
+ else
+ {
+ return LW_FALSE;
+ }
}
return LW_FALSE;
}
@@ -3434,6 +3467,7 @@
/* Are A-plane and B-plane basically the same? */
ab_dot = dot_product(&AN, &BN);
+
if ( FP_EQUALS(fabs(ab_dot), 1.0) )
{
/* Co-linear case */
Modified: trunk/liblwgeom/lwgeodetic_tree.c
===================================================================
--- trunk/liblwgeom/lwgeodetic_tree.c 2019-01-11 13:50:49 UTC (rev 17136)
+++ trunk/liblwgeom/lwgeodetic_tree.c 2019-01-11 21:07:17 UTC (rev 17137)
@@ -479,7 +479,7 @@
* KNOWN PROBLEM: Grazings (think of a sharp point, just touching the
* stabline) will be counted for one, which will throw off the count.
*/
-int circ_tree_contains_point(const CIRC_NODE* node, const POINT2D* pt, const POINT2D* pt_outside, int* on_boundary)
+int circ_tree_contains_point(const CIRC_NODE* node, const POINT2D* pt, const POINT2D* pt_outside, int level, int* on_boundary)
{
GEOGRAPHIC_POINT closest;
GEOGRAPHIC_EDGE stab_edge, edge;
@@ -493,7 +493,7 @@
geog2cart(&(stab_edge.start), &S1);
geog2cart(&(stab_edge.end), &S2);
- LWDEBUG(3, "entered");
+ LWDEBUGF(3, "%*s entered", level, "");
/*
* If the stabline doesn't cross within the radius of a node, there's no
@@ -500,18 +500,18 @@
* way it can cross.
*/
- LWDEBUGF(3, "working on node %p, edge_num %d, radius %g, center POINT(%g %g)", node, node->edge_num, node->radius, rad2deg(node->center.lon), rad2deg(node->center.lat));
+ LWDEBUGF(3, "%*s :working on node %p, edge_num %d, radius %g, center POINT(%.12g %.12g)", level, "", node, node->edge_num, node->radius, rad2deg(node->center.lon), rad2deg(node->center.lat));
d = edge_distance_to_point(&stab_edge, &(node->center), &closest);
- LWDEBUGF(3, "edge_distance_to_point=%g, node_radius=%g", d, node->radius);
+ LWDEBUGF(3, "%*s :edge_distance_to_point=%g, node_radius=%g", level, "", d, node->radius);
if ( FP_LTEQ(d, node->radius) )
{
- LWDEBUGF(3,"entering this branch (%p)", node);
+ LWDEBUGF(3,"%*s :entering this branch (%p)", level, "", node);
/* Return the crossing number of this leaf */
if ( circ_node_is_leaf(node) )
{
int inter;
- LWDEBUGF(3, "leaf node calculation (edge %d)", node->edge_num);
+ LWDEBUGF(3, "%*s :leaf node calculation (edge %d)", level, "", node->edge_num);
geographic_point_init(node->p1->x, node->p1->y, &(edge.start));
geographic_point_init(node->p2->x, node->p2->y, &(edge.end));
geog2cart(&(edge.start), &E1);
@@ -518,24 +518,41 @@
geog2cart(&(edge.end), &E2);
inter = edge_intersects(&S1, &S2, &E1, &E2);
+ LWDEBUGF(3, "%*s :inter = %d", level, "", inter);
if ( inter & PIR_INTERSECTS )
{
- LWDEBUG(3," got stab line edge_intersection with this edge!");
+ LWDEBUGF(3,"%*s ::got stab line edge_intersection with this edge!", level, "");
/* To avoid double counting crossings-at-a-vertex, */
/* always ignore crossings at "lower" ends of edges*/
+ GEOGRAPHIC_POINT e1, e2;
+ cart2geog(&E1,&e1); cart2geog(&E2,&e2);
+ LWDEBUGF(3,"%*s LINESTRING(%.15g %.15g,%.15g %.15g)", level, "",
+ pt->x, pt->y,
+ pt_outside->x, pt_outside->y
+ );
+
+ LWDEBUGF(3,"%*s LINESTRING(%.15g %.15g,%.15g %.15g)", level, "",
+ rad2deg(e1.lon), rad2deg(e1.lat),
+ rad2deg(e2.lon), rad2deg(e2.lat)
+ );
+
if ( inter & PIR_B_TOUCH_RIGHT || inter & PIR_COLINEAR )
{
- LWDEBUG(3," rejecting stab line grazing by left-side edge");
+ LWDEBUGF(3,"%*s ::rejecting stab line grazing by left-side edge", level, "");
return 0;
}
else
{
- LWDEBUG(3," accepting stab line intersection");
+ LWDEBUGF(3,"%*s ::accepting stab line intersection", level, "");
return 1;
}
}
+ else
+ {
+ LWDEBUGF(3,"%*s edge does not intersect", level, "");
+ }
}
/* Or, add up the crossing numbers of all children of this node. */
else
@@ -543,9 +560,8 @@
c = 0;
for ( i = 0; i < node->num_nodes; i++ )
{
- LWDEBUG(3,"internal node calculation");
- LWDEBUGF(3," calling circ_tree_contains_point on child %d!", i);
- c += circ_tree_contains_point(node->nodes[i], pt, pt_outside, on_boundary);
+ LWDEBUGF(3,"%*s calling circ_tree_contains_point on child %d!", level, "", i);
+ c += circ_tree_contains_point(node->nodes[i], pt, pt_outside, level + 1, on_boundary);
}
return c % 2;
}
@@ -552,9 +568,8 @@
}
else
{
- LWDEBUGF(3,"skipping this branch (%p)", node);
+ LWDEBUGF(3,"%*s skipping this branch (%p)", level, "", node);
}
-
return 0;
}
@@ -654,11 +669,12 @@
uint32_t i;
LWDEBUGF(4, "entered, min_dist=%.8g max_dist=%.8g, type1=%d, type2=%d", *min_dist, *max_dist, n1->geom_type, n2->geom_type);
-/*
- circ_tree_print(n1, 0);
- circ_tree_print(n2, 0);
-*/
+ // printf("-==-\n");
+ // circ_tree_print(n1, 0);
+ // printf("--\n");
+ // circ_tree_print(n2, 0);
+
/* Short circuit if we've already hit the minimum */
if( *min_dist < threshold || *min_dist == 0.0 )
return *min_dist;
@@ -683,7 +699,7 @@
POINT2D pt;
circ_tree_get_point(n2, &pt);
LWDEBUGF(4, "n1 is polygon, testing if contains (%.5g,%.5g)", pt.x, pt.y);
- if ( circ_tree_contains_point(n1, &pt, &(n1->pt_outside), NULL) )
+ if ( circ_tree_contains_point(n1, &pt, &(n1->pt_outside), 0, NULL) )
{
LWDEBUG(4, "it does");
*min_dist = 0.0;
@@ -699,7 +715,7 @@
POINT2D pt;
circ_tree_get_point(n1, &pt);
LWDEBUGF(4, "n2 is polygon, testing if contains (%.5g,%.5g)", pt.x, pt.y);
- if ( circ_tree_contains_point(n2, &pt, &(n2->pt_outside), NULL) )
+ if ( circ_tree_contains_point(n2, &pt, &(n2->pt_outside), 0, NULL) )
{
LWDEBUG(4, "it does");
geographic_point_init(pt.x, pt.y, closest1);
@@ -791,6 +807,7 @@
/* tests above. */
if ( n1->geom_type && lwtype_is_collection(n1->geom_type) )
{
+ printf("arg1 is a collection\n");
circ_internal_nodes_sort(n1->nodes, n1->num_nodes, n2);
for ( i = 0; i < n1->num_nodes; i++ )
{
@@ -800,6 +817,7 @@
}
else if ( n2->geom_type && lwtype_is_collection(n2->geom_type) )
{
+ printf("arg2 is a collection\n");
circ_internal_nodes_sort(n2->nodes, n2->num_nodes, n1);
for ( i = 0; i < n2->num_nodes; i++ )
{
@@ -875,7 +893,7 @@
}
if ( node->geom_type == POLYGONTYPE )
{
- printf(" O(%.5g %.5g)", node->pt_outside.x, node->pt_outside.y);
+ printf(" O(%.15g %.15g)", node->pt_outside.x, node->pt_outside.y);
}
printf("\n");
}
Modified: trunk/liblwgeom/lwgeodetic_tree.h
===================================================================
--- trunk/liblwgeom/lwgeodetic_tree.h 2019-01-11 13:50:49 UTC (rev 17136)
+++ trunk/liblwgeom/lwgeodetic_tree.h 2019-01-11 21:07:17 UTC (rev 17137)
@@ -50,7 +50,7 @@
void circ_tree_print(const CIRC_NODE* node, int depth);
CIRC_NODE* circ_tree_new(const POINTARRAY* pa);
void circ_tree_free(CIRC_NODE* node);
-int circ_tree_contains_point(const CIRC_NODE* node, const POINT2D* pt, const POINT2D* pt_outside, int* on_boundary);
+int circ_tree_contains_point(const CIRC_NODE* node, const POINT2D* pt, const POINT2D* pt_outside, int level, int* on_boundary);
double circ_tree_distance_tree(const CIRC_NODE* n1, const CIRC_NODE* n2, const SPHEROID *spheroid, double threshold);
CIRC_NODE* lwgeom_calculate_circ_tree(const LWGEOM* lwgeom);
int circ_tree_get_point(const CIRC_NODE* node, POINT2D* pt);
More information about the postgis-tickets
mailing list