Mercurial > gemma
changeset 5073:9bbe9d0a2dfb time-sliding
merge with default
author | Thomas Junk <thomas.junk@intevation.de> |
---|---|
date | Thu, 12 Mar 2020 08:19:49 +0100 |
parents | 51b30f7c64a2 (current diff) ae3a1392f9d0 (diff) |
children | 8ec4f1a0db86 |
files | |
diffstat | 11 files changed, 258 insertions(+), 98 deletions(-) [+] |
line wrap: on
line diff
--- a/pkg/imports/fm.go Wed Mar 11 16:45:13 2020 +0100 +++ b/pkg/imports/fm.go Thu Mar 12 08:19:49 2020 +0100 @@ -219,7 +219,7 @@ ), consume, createInvalidation("bcnlat_hydro"), - newPointSlice(func() interface{} { return new(bcnlatHydroProperties) }), + newPointFeature(func() interface{} { return new(bcnlatHydroProperties) }), ), }) @@ -235,7 +235,7 @@ ), consume, createInvalidation("bcnlat_ienc"), - newPointSlice(func() interface{} { return new(bcnlatIencProperties) }), + newPointFeature(func() interface{} { return new(bcnlatIencProperties) }), ), }) @@ -251,7 +251,7 @@ ), consume, createInvalidation("boylat_hydro"), - newPointSlice(func() interface{} { return new(boylatHydroProperties) }), + newPointFeature(func() interface{} { return new(boylatHydroProperties) }), ), }) @@ -267,7 +267,7 @@ ), consume, createInvalidation("boylat_ienc"), - newPointSlice(func() interface{} { return new(boylatIencProperties) }), + newPointFeature(func() interface{} { return new(boylatIencProperties) }), ), }) @@ -283,7 +283,7 @@ ), consume, createInvalidation("boycar"), - newPointSlice(func() interface{} { return new(boycarProperties) }), + newPointFeature(func() interface{} { return new(boycarProperties) }), ), }) @@ -298,7 +298,7 @@ ), consume, createInvalidation("boysaw"), - newPointSlice(func() interface{} { return new(boysawProperties) }), + newPointFeature(func() interface{} { return new(boysawProperties) }), ), }) @@ -314,7 +314,7 @@ ), consume, createInvalidation("boyspp"), - newPointSlice(func() interface{} { return new(boysppProperties) }), + newPointFeature(func() interface{} { return new(boysppProperties) }), ), }) @@ -329,7 +329,7 @@ ), consume, createInvalidation("daymar_hydro"), - newPointSlice(func() interface{} { return new(daymarHydroProperties) }), + newPointFeature(func() interface{} { return new(daymarHydroProperties) }), ), }) @@ -345,7 +345,7 @@ ), consume, createInvalidation("daymar_ienc"), - newPointSlice(func() interface{} { return new(daymarIencProperties) }), + newPointFeature(func() interface{} { return new(daymarIencProperties) }), ), }) @@ -364,7 +364,7 @@ ), consume, createInvalidation("lights"), - newPointSlice(func() interface{} { return new(lightsProperties) }), + newPointFeature(func() interface{} { return new(lightsProperties) }), ), }) @@ -383,7 +383,7 @@ ), consume, createInvalidation("notmrk"), - newPointSlice(func() interface{} { return new(notmrkProperties) }), + newPointFeature(func() interface{} { return new(notmrkProperties) }), ), }) @@ -398,7 +398,7 @@ ), consume, createInvalidation("rtpbcn"), - newPointSlice(func() interface{} { return new(rtpbcnProperties) }), + newPointFeature(func() interface{} { return new(rtpbcnProperties) }), ), }) @@ -413,7 +413,7 @@ ), consume, createInvalidation("topmar"), - newPointSlice(func() interface{} { return new(topmarProperties) }), + newPointFeature(func() interface{} { return new(topmarProperties) }), ), }) } @@ -567,7 +567,7 @@ func consume( spc *SQLGeometryConsumer, - points, properties interface{}, + geom, properties interface{}, epsg int, ) error { var fmid int64 @@ -576,7 +576,7 @@ spc.ctx, append( []interface{}{ - points.(*pointSlice).asWKB(), + geom.(interface{ asWKB() []byte }).asWKB(), epsg, }, structs.Values(properties)...)...,
--- a/pkg/imports/wfsjob.go Wed Mar 11 16:45:13 2020 +0100 +++ b/pkg/imports/wfsjob.go Thu Mar 12 08:19:49 2020 +0100 @@ -38,7 +38,7 @@ Commit() error Rollback() error - NewFeature() (kind string, geom interface{}, properties interface{}) + NewFeature() (kind string, properties interface{}) Consume(geom, properties interface{}, epsg int) error } @@ -58,6 +58,22 @@ } ) +var ( + kindToGeometry = map[string]func() interface{}{ + // TODO: extend me! + "Point": func() interface{} { return new(pointSlice) }, + "LineString": func() interface{} { return new(lineSlice) }, + "MultiLineString": func() interface{} { return new(multiLineSlice) }, + } + + wrapGeomKind = map[[2]string]func(interface{}) interface{}{ + // TODO: extend me! + {"LineString", "MultiLineString"}: func(x interface{}) interface{} { + return &multiLineSlice{*x.(*lineSlice)} + }, + } +) + func (wfjc *WFSFeatureJobCreator) Description() string { return wfjc.description } @@ -176,31 +192,44 @@ continue } - kind, geom, props := consumer.NewFeature() + kind, props := consumer.NewFeature() if err := json.Unmarshal(*feature.Properties, props); err != nil { badProperties++ continue } - if feature.Geometry.Type == kind { - if err := json.Unmarshal(*feature.Geometry.Coordinates, geom); err != nil { - return err - } + // Look if we can deserialize given type + makeGeom := kindToGeometry[feature.Geometry.Type] + if makeGeom == nil { + unsupported[feature.Geometry.Type]++ + continue + } - err := consumer.Consume(geom, props, epsg) - switch { - case err == ErrFeatureDuplicated: - dupes++ - case err == ErrFeatureIgnored: - // be silent - case err != nil: - return err - default: - features++ + // Optional wrapping + wrap := func(x interface{}) interface{} { return x } + if feature.Geometry.Type != kind { + // Look if we can wrap it + if wrap = wrapGeomKind[[2]string{feature.Geometry.Type, kind}]; wrap == nil { + unsupported[feature.Geometry.Type]++ + continue } - } else { - unsupported[feature.Geometry.Type]++ + } + + geom := makeGeom() + if err := json.Unmarshal(*feature.Geometry.Coordinates, geom); err != nil { + return err + } + + switch err := consumer.Consume(wrap(geom), props, epsg); { + case err == ErrFeatureDuplicated: + dupes++ + case err == ErrFeatureIgnored: + // be silent + case err != nil: + return err + default: + features++ } } return nil @@ -251,7 +280,7 @@ tx *sql.Tx feedback Feedback consume func(*SQLGeometryConsumer, interface{}, interface{}, int) error - newFeature func() (string, interface{}, interface{}) + newFeature func() (string, interface{}) preCommit func(*SQLGeometryConsumer) error savepoint func(func() error) error stmts []*sql.Stmt @@ -285,7 +314,7 @@ return err } -func (sgc *SQLGeometryConsumer) NewFeature() (string, interface{}, interface{}) { +func (sgc *SQLGeometryConsumer) NewFeature() (string, interface{}) { return sgc.newFeature() } @@ -308,7 +337,7 @@ init func(*SQLGeometryConsumer) error, consume func(*SQLGeometryConsumer, interface{}, interface{}, int) error, preCommit func(*SQLGeometryConsumer) error, - newFeature func() (string, interface{}, interface{}), + newFeature func() (string, interface{}), ) func(context.Context, *sql.Conn, Feedback) (WFSFeatureConsumer, error) { return func(ctx context.Context, conn *sql.Conn, feedback Feedback) (WFSFeatureConsumer, error) {
--- a/pkg/imports/wkb.go Wed Mar 11 16:45:13 2020 +0100 +++ b/pkg/imports/wkb.go Thu Mar 12 08:19:49 2020 +0100 @@ -25,23 +25,17 @@ ) type ( - pointSlice []float64 - lineSlice [][]float64 - polygonSlice [][][]float64 + pointSlice []float64 + lineSlice [][]float64 + multiLineSlice []lineSlice + polygonSlice [][][]float64 ) -func newPointSlice(newProperties func() interface{}) func() (string, interface{}, interface{}) { - return func() (string, interface{}, interface{}) { - return "Point", new(pointSlice), newProperties() - } +func newPointFeature(newProperties func() interface{}) func() (string, interface{}) { + return func() (string, interface{}) { return "Point", newProperties() } } -func (ls lineSlice) asWKB() []byte { - - size := 1 + 4 + 4 + len(ls)*(2*8) - - buf := bytes.NewBuffer(make([]byte, 0, size)) - +func (ls lineSlice) toWKB(buf *bytes.Buffer) { binary.Write(buf, binary.LittleEndian, wkb.NDR) binary.Write(buf, binary.LittleEndian, wkb.LineString) binary.Write(buf, binary.LittleEndian, uint32(len(ls))) @@ -57,6 +51,14 @@ binary.Write(buf, binary.LittleEndian, math.Float64bits(lat)) binary.Write(buf, binary.LittleEndian, math.Float64bits(lon)) } +} + +func (ls lineSlice) asWKB() []byte { + + size := 1 + 4 + 4 + len(ls)*(2*8) + + buf := bytes.NewBuffer(make([]byte, 0, size)) + ls.toWKB(buf) return buf.Bytes() } @@ -70,6 +72,27 @@ return lr } +func (mls multiLineSlice) asWKB() []byte { + + size := 1 + 4 + 4 + + for _, ls := range mls { + size += 1 + 4 + 4 + len(ls)*(2*8) + } + + buf := bytes.NewBuffer(make([]byte, 0, size)) + + binary.Write(buf, binary.LittleEndian, wkb.NDR) + binary.Write(buf, binary.LittleEndian, wkb.MultiLineString) + binary.Write(buf, binary.LittleEndian, uint32(len(mls))) + + for _, ls := range mls { + ls.toWKB(buf) + } + + return buf.Bytes() +} + func (p pointSlice) asWKB() []byte { size := 1 + 4 + 2*8
--- a/pkg/imports/wx.go Wed Mar 11 16:45:13 2020 +0100 +++ b/pkg/imports/wx.go Thu Mar 12 08:19:49 2020 +0100 @@ -100,16 +100,18 @@ SELECT users.current_user_area_utm() AS a ) INSERT INTO waterway.waterway_axis (wtwaxs, objnam, nobjnam) -SELECT dmp.geom, $3, $4 +SELECT + ST_Multi(ST_Node(ST_CollectionExtract(ST_Transform(new_ax, 4326), 2))), + $3, $4 FROM ST_GeomFromWKB($1, $2::integer) AS new_line (new_line), - ST_Dump(ST_Transform(ST_CollectionExtract( + LATERAL (SELECT CASE WHEN pg_has_role('sys_admin', 'MEMBER') - THEN ST_Node(ST_Transform(new_line, - best_utm(ST_Transform(new_line, 4326)))) + THEN new_line ELSE ST_Intersection((SELECT a FROM resp), ST_Node(ST_Transform(new_line, (SELECT ST_SRID(a) FROM resp)))) - END, - 2), 4326)) AS dmp + END) AS new_ax (new_ax) + -- Do nothing if intersection is empty: + WHERE NOT ST_IsEmpty(new_ax) RETURNING id ` ) @@ -222,47 +224,34 @@ nobjnam = sql.NullString{String: *props.NObjNnm, Valid: true} } + var ls multiLineSlice switch feature.Geometry.Type { case "LineString": var l lineSlice if err := json.Unmarshal(*feature.Geometry.Coordinates, &l); err != nil { return err } - if err := storeLinestring( - ctx, - savepoint, - feedback, - l, - epsg, - props, - nobjnam, - &outside, - &features, - insertStmt); err != nil { - return err - } + ls = append(ls, l) case "MultiLineString": - var ls []lineSlice if err := json.Unmarshal(*feature.Geometry.Coordinates, &ls); err != nil { return err } - for _, l := range ls { - if err := storeLinestring( - ctx, - savepoint, - feedback, - l, - epsg, - props, - nobjnam, - &outside, - &features, - insertStmt); err != nil { - return err - } - } default: unsupported[feature.Geometry.Type]++ + continue + } + if err := storeLinestring( + ctx, + savepoint, + feedback, + ls, + epsg, + props, + nobjnam, + &outside, + &features, + insertStmt); err != nil { + return err } } return nil @@ -302,7 +291,7 @@ ctx context.Context, savepoint func(func() error) error, feedback Feedback, - l lineSlice, + l multiLineSlice, epsg int, props waterwayAxisProperties, nobjnam sql.NullString,
--- a/schema/gemma.sql Wed Mar 11 16:45:13 2020 +0100 +++ b/schema/gemma.sql Thu Mar 12 08:19:49 2020 +0100 @@ -620,7 +620,7 @@ CREATE TABLE waterway_axis ( id int PRIMARY KEY GENERATED BY DEFAULT AS IDENTITY, - wtwaxs geography(LINESTRING, 4326) NOT NULL + wtwaxs geography(MULTILINESTRING, 4326) NOT NULL CHECK(ST_IsSimple(CAST(wtwaxs AS geometry))), -- TODO: Do we need to check data set quality (DRC 2.1.6)? objnam varchar NOT NULL,
--- a/schema/gemma_tests.sql Wed Mar 11 16:45:13 2020 +0100 +++ b/schema/gemma_tests.sql Thu Mar 12 08:19:49 2020 +0100 @@ -26,8 +26,8 @@ SELECT throws_ok($$ INSERT INTO waterway.waterway_axis (wtwaxs, objnam) VALUES - (ST_GeogFromText('LINESTRING(0 0, 1 1)'), 'test'), - (ST_GeogFromText('LINESTRING(0 0, 1 1)'), 'test') + (ST_GeogFromText('MULTILINESTRING((0 0, 1 1))'), 'test'), + (ST_GeogFromText('MULTILINESTRING((0 0, 1 1))'), 'test') $$, 23505, NULL, 'No duplicate geometries can be inserted into waterway_axis');
--- a/schema/isrs_functions.sql Wed Mar 11 16:45:13 2020 +0100 +++ b/schema/isrs_functions.sql Thu Mar 12 08:19:49 2020 +0100 @@ -61,9 +61,11 @@ z = best_utm(stretch); CREATE TEMP TABLE axis AS - SELECT id, wtwaxs, ST_Boundary(wtwaxs) AS bdr - FROM (SELECT id, ST_Transform(wtwaxs::geometry, z) AS wtwaxs - FROM waterway.waterway_axis) AS axs; + SELECT row_number() OVER () AS id, + geom AS wtwaxs, + ST_Boundary(geom) AS bdr + FROM waterway.waterway_axis, + ST_Dump(ST_Transform(wtwaxs::geometry, z)); CREATE INDEX axs_bdr ON axis USING GiST (bdr); ANALYZE axis;
--- a/schema/tap_tests_data.sql Wed Mar 11 16:45:13 2020 +0100 +++ b/schema/tap_tests_data.sql Thu Mar 12 08:19:49 2020 +0100 @@ -124,22 +124,23 @@ ); INSERT INTO waterway.waterway_axis (wtwaxs, objnam) VALUES ( - ST_SetSRID(ST_CurveToLine( - 'CIRCULARSTRING(0 0, 0.5 0.5, 0.6 0.4)'), + ST_SetSRID( + ST_Multi(ST_CurveToLine('CIRCULARSTRING(0 0, 0.5 0.5, 0.6 0.4)')), 4326), 'testriver' ), ( - ST_SetSRID(ST_CurveToLine('CIRCULARSTRING(0.6 0.4, 1 0, 1.5 -0.00001)'), + ST_SetSRID( + ST_Multi(ST_CurveToLine('CIRCULARSTRING(0.6 0.4, 1 0, 1.5 -0.00001)')), 4326), 'testriver' ), ( - ST_SetSRID('LINESTRING(0.5 0.5, 1 1)'::geometry, 4326), + ST_SetSRID('MULTILINESTRING((0.5 0.5, 1 1))'::geometry, 4326), 'testriver' ), ( - ST_SetSRID('LINESTRING(1.5 0, 1.55001 0)'::geometry, 4326), + ST_SetSRID('MULTILINESTRING((1.5 0, 1.55001 0))'::geometry, 4326), 'testriver' ), ( - ST_SetSRID('LINESTRING(1.55 0, 2 0)'::geometry, 4326), + ST_SetSRID('MULTILINESTRING((1.55 0, 2 0))'::geometry, 4326), 'testriver' );
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/schema/updates/1423/01.axis_as_multilinestring.sql Thu Mar 12 08:19:49 2020 +0100 @@ -0,0 +1,11 @@ +-- Cannot alter type of a column used in a trigger definition +DROP TRIGGER waterway_axis_wtwaxs_unique ON waterway.waterway_axis; + +ALTER TABLE waterway.waterway_axis + ALTER wtwaxs TYPE geography(MULTILINESTRING, 4326) + USING ST_Multi(CAST(wtwaxs AS geometry)); + +-- Re-create trigger +CREATE CONSTRAINT TRIGGER waterway_axis_wtwaxs_unique + AFTER INSERT OR UPDATE OF wtwaxs ON waterway.waterway_axis + FOR EACH ROW EXECUTE FUNCTION prevent_st_equals('wtwaxs');
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/schema/updates/1423/02.adapt_func.sql Thu Mar 12 08:19:49 2020 +0100 @@ -0,0 +1,105 @@ +CREATE OR REPLACE FUNCTION ISRSrange_axis( + stretch isrsrange, + tolerance float + -- in m, up to which linestrings will be connected at their boundary +) RETURNS geometry +AS $$ +DECLARE z int; +DECLARE result_geom geometry; +BEGIN + -- Find best matchting UTM zone + z = best_utm(stretch); + + CREATE TEMP TABLE axis AS + SELECT row_number() OVER () AS id, + geom AS wtwaxs, + ST_Boundary(geom) AS bdr + FROM waterway.waterway_axis, + ST_Dump(ST_Transform(wtwaxs::geometry, z)); + CREATE INDEX axs_bdr ON axis USING GiST (bdr); + ANALYZE axis; + + WITH RECURSIVE + -- In order to guarantee the following ST_Covers to work, + -- snap distance mark coordinates to axis + points0 AS ( + SELECT ST_ClosestPoint( + wtwaxs, + ST_Transform(geom, z)) AS geom + FROM ST_Dump(ISRSrange_points(stretch)), ( + SELECT ST_Collect(wtwaxs) AS wtwaxs + FROM axis) AS ax), + -- Ensure two distinct points on axis have been found + points AS ( + SELECT geom + FROM points0 + WHERE 2 = (SELECT count(DISTINCT geom) FROM points0)), + axis_snapped AS ( + -- Iteratively connect non-contiguous axis chunks + -- to find the contiguous axis on which given distance marks lie + (SELECT ARRAY[id] AS ids, wtwaxs + FROM axis, points + WHERE ST_Intersects( + ST_Buffer(axis.wtwaxs, 0.0001), points.geom) + FETCH FIRST ROW ONLY) + UNION + -- Connect endpoint of next linestring with closest + -- endpoint of merged linestring until a contiguous + -- linestring connecting both distance marks is build up + (SELECT refids || id, + ST_LineMerge(ST_Collect(ARRAY( + -- Linestring build up so far + SELECT refgeom + UNION + -- Fill eventual gap + SELECT ST_ShortestLine( + ST_Boundary(refgeom), bdr) + UNION + -- Linestring to be added + SELECT geom))) + FROM axis_snapped AS axis_snapped (refids, refgeom), + axis AS axis (id, geom, bdr), + (SELECT ST_Collect(points.geom) AS pts + FROM points) AS points + WHERE id <> ALL(refids) + AND ST_DWithin( + ST_Boundary(refgeom), bdr, tolerance) + AND NOT ST_Covers(ST_Buffer(refgeom, 0.0001), points.pts) + ORDER BY ST_Boundary(refgeom) <-> bdr + FETCH FIRST ROW ONLY)), + axis_segment AS ( + -- Fetch end result from snapping + SELECT wtwaxs AS line + FROM axis_snapped, + (SELECT ST_Collect(points.geom) AS pts + FROM points) AS points + -- Return end result only if both distance marks were connected + WHERE ST_Covers(ST_Buffer(wtwaxs, 0.0001), points.pts)) + -- Use linear referencing to clip axis between distance marks. + -- Simplification is used to work-around the problem, that + -- ST_LineSubstring might generate very small line segments at an + -- end of the resulting linestring, that significantly differ from + -- the direction of the input linestring due to finite precision + -- of the calculation. The generated small segment of the + -- resulting line would lead e.g. to unexpected results in an area + -- generated by ISRSrange_area(). + SELECT ST_SimplifyPreserveTopology(ST_LineSubstring( + axis_segment.line, min(fractions.f), max(fractions.f)), + 0.0001) AS line + INTO STRICT result_geom + FROM axis_segment, LATERAL ( + SELECT ST_LineLocatePoint(axis_segment.line, points.geom) AS f + FROM points) AS fractions + GROUP BY axis_segment.line; + + -- Drop temporary table to avoid side effects on PostgreSQL's MVCC, + -- because otherwise subsequent invocations of the function will not see + -- changes on the underlying waterway.waterway_axis that might have + -- occured. + DROP TABLE axis; + + RETURN result_geom; +END; + $$ + LANGUAGE plpgsql + PARALLEL RESTRICTED;