Mercurial > gemma
view pkg/models/cross.go @ 721:ca82698349b7
Merged octree branch into default.
author | Sascha L. Teichmann <sascha.teichmann@intevation.de> |
---|---|
date | Fri, 21 Sep 2018 17:37:34 +0200 |
parents | bd215c4325ce |
children | ce23a4192788 |
line wrap: on
line source
package models import ( "bytes" "encoding/binary" "encoding/json" "errors" "fmt" "math" "sort" "time" ) type ( GeoJSONFeature string GeoJSONLineStringType string GeoJSONDate struct{ time.Time } GeoJSONLineString struct { Type GeoJSONLineStringType `json:"type"` } GeoJSONCoordinate struct { Lat float64 Lon float64 } GeoJSONCoordinateZ struct { Lat float64 Lon float64 Z float64 } CrossSectionInputProperties struct { Bottleneck string `json:"bottleneck"` Date GeoJSONDate `json:"date"` } GeoJSONLineCoordinates []GeoJSONCoordinate GeoJSONLineCoordinatesZ []GeoJSONCoordinateZ GeoJSONMultiLineCoordinatesZ []GeoJSONLineCoordinatesZ CrossSectionInputGeometry struct { Type GeoJSONLineStringType `json:"type"` Coordinates GeoJSONLineCoordinates `json:"coordinates"` } CrossSectionInput struct { Type GeoJSONFeature `json:"type"` Geometry CrossSectionInputGeometry `json:"geometry"` Properties CrossSectionInputProperties `json:"properties"` } CrossSectionOutputGeometry struct { Type string `json:"type"` Coordinates GeoJSONMultiLineCoordinatesZ `json:"coordinates"` } CrossSectionOutput struct { Type string `json:"type"` Geometry CrossSectionOutputGeometry `json:"geometry"` Properties map[string]interface{} `json:"properties"` } ) var ( errNoGeoJSONFeature = errors.New("Not a valid GeoJSON feature") errNoGeoJSONLineStringType = errors.New("Not a valid GeoJSON linestring type") errTooLessComponents = errors.New("Too less components in coordinate") errTooLessCoordinates = errors.New("Too less coordinates") ) func (lc *GeoJSONLineCoordinates) UnmarshalJSON(data []byte) error { var coords []GeoJSONCoordinate if err := json.Unmarshal(data, &coords); err != nil { return err } if len(coords) < 2 { return errTooLessCoordinates } *lc = GeoJSONLineCoordinates(coords) return nil } func (d *GeoJSONDate) UnmarshalJSON(data []byte) error { var s string if err := json.Unmarshal(data, &s); err != nil { return err } t, err := time.Parse("2006-01-02", s) if err != nil { return err } *d = GeoJSONDate{t} return nil } func (c *GeoJSONCoordinate) UnmarshalJSON(data []byte) error { var pos []float64 if err := json.Unmarshal(data, &pos); err != nil { return err } if len(pos) < 2 { return errTooLessComponents } *c = GeoJSONCoordinate{Lat: pos[0], Lon: pos[1]} return nil } func (t *GeoJSONFeature) UnmarshalJSON(data []byte) error { var s string if err := json.Unmarshal(data, &s); err != nil { return err } if s != "Feature" { return errNoGeoJSONFeature } *t = GeoJSONFeature(s) return nil } func (t *GeoJSONLineStringType) UnmarshalJSON(data []byte) error { var s string if err := json.Unmarshal(data, &s); err != nil { return err } if s != "LineString" { return errNoGeoJSONLineStringType } *t = GeoJSONLineStringType(s) return nil } const ( wkbXDR byte = 0 wkbNDR byte = 1 wkbPoint uint32 = 1 wkbPointZ uint32 = 1000 + 1 wkbLineString uint32 = 2 wkbLineStringZ uint32 = 1000 + 2 ) func (lc GeoJSONLineCoordinates) AsWKB() []byte { size := 1 + 4 + 4 + len(lc)*(1+4+2*8) buf := bytes.NewBuffer(make([]byte, 0, size)) binary.Write(buf, binary.LittleEndian, wkbNDR) binary.Write(buf, binary.LittleEndian, wkbLineString) binary.Write(buf, binary.LittleEndian, uint32(len(lc))) for i := range lc { c := &lc[i] binary.Write(buf, binary.LittleEndian, math.Float64bits(c.Lat)) binary.Write(buf, binary.LittleEndian, math.Float64bits(c.Lon)) } return buf.Bytes() } func (cz GeoJSONCoordinateZ) MarshalJSON() ([]byte, error) { var buf bytes.Buffer fmt.Fprintf(&buf, "[%.8f,%.8f,%.8f]", cz.Lat, cz.Lon, cz.Z) return buf.Bytes(), nil } func (lcz *GeoJSONLineCoordinatesZ) Scan(src interface{}) error { data, ok := src.([]byte) if !ok { return errNoByteSlice } return lcz.FromWKB(data) } func (lcz *GeoJSONLineCoordinatesZ) FromWKB(data []byte) error { r := bytes.NewReader(data) endian, err := r.ReadByte() var order binary.ByteOrder switch { case err != nil: return err case endian == wkbNDR: order = binary.LittleEndian case endian == wkbXDR: order = binary.BigEndian default: return fmt.Errorf("unknown byte order %x", endian) } var geomType uint32 err = binary.Read(r, order, &geomType) switch { case err != nil: return err case geomType != wkbLineStringZ: return fmt.Errorf("unknown geometry type %x", geomType) } var num uint32 if err = binary.Read(r, order, &num); err != nil { return err } coords := make(GeoJSONLineCoordinatesZ, num) for i := range coords { c := &coords[i] for _, addr := range []*float64{&c.Lat, &c.Lon, &c.Z} { if err = binary.Read(r, order, addr); err != nil { return err } } } *lcz = coords return nil } func (cz GeoJSONCoordinateZ) Equals(other GeoJSONCoordinateZ) bool { const ( xyEps = 1e-7 zEps = 1e-5 ) return math.Abs(cz.Lat-other.Lat) < xyEps && math.Abs(cz.Lon-other.Lon) < xyEps && math.Abs(cz.Z-other.Z) < zEps } func (cz GeoJSONCoordinateZ) mid(other GeoJSONCoordinateZ) GeoJSONCoordinateZ { return GeoJSONCoordinateZ{ Lon: (cz.Lon + other.Lon) * 0.5, Lat: (cz.Lat + other.Lat) * 0.5, Z: (cz.Z + other.Z) * 0.5, } } func deg2rad(d float64) float64 { return d * math.Pi / 180.0 } func (cz GeoJSONCoordinateZ) Distance(other GeoJSONCoordinateZ) float64 { const EarthRadius = 6378137.0 dLat := deg2rad(other.Lat - cz.Lat) dLng := math.Abs(deg2rad(other.Lon - cz.Lon)) if dLng > math.Pi { dLng = 2*math.Pi - dLng } x := dLng * math.Cos(deg2rad((cz.Lat+other.Lat)/2.0)) return math.Sqrt(dLat*dLat+x*x) * EarthRadius } type point struct { x float64 y float64 } func (cz GeoJSONCoordinateZ) quant() point { const xyScale = 1e6 return point{ math.Round(cz.Lon*xyScale) / xyScale, math.Round(cz.Lat*xyScale) / xyScale, } } func (lcz GeoJSONLineCoordinatesZ) join(other GeoJSONLineCoordinatesZ) GeoJSONLineCoordinatesZ { l1, l2 := len(lcz), len(other) n := make(GeoJSONLineCoordinatesZ, l1+l2-1) copy(n, lcz[:l1-1]) n[l1-1] = lcz[l1-1].mid(other[0]) copy(n[l1:], other[1:]) return n } func (ml GeoJSONMultiLineCoordinatesZ) Join() GeoJSONMultiLineCoordinatesZ { if len(ml) < 2 { return ml } type value struct { coords GeoJSONLineCoordinatesZ order int } heads := make(map[point]*value) for i, coords := range ml { key := coords[len(coords)-1].quant() if v := heads[key]; v != nil { delete(heads, key) v.coords = coords.join(v.coords) heads[v.coords[0].quant()] = v } else { heads[coords[0].quant()] = &value{coords, i} } } again: for k, v1 := range heads { key := v1.coords[len(v1.coords)-1].quant() if k == key { // cycle detected continue } if v2 := heads[key]; v2 != nil { delete(heads, key) v1.coords = v1.coords.join(v2.coords) goto again } } // To make it deterministic sort in input order. tmp := make([]*value, len(heads)) var i int for _, v := range heads { tmp[i] = v i++ } sort.Slice(tmp, func(i, j int) bool { return tmp[i].order < tmp[j].order }) out := make(GeoJSONMultiLineCoordinatesZ, len(tmp)) for i, v := range tmp { out[i] = v.coords } return out } /* func (ml GeoJSONMultiLineCoordinatesZ) Join() GeoJSONMultiLineCoordinatesZ { if len(ml) < 2 { return ml } lists := make(GeoJSONMultiLineCoordinatesZ, len(ml)) copy(lists, ml) next: for j := 0; j < len(lists); { front := lists[j] tail := front[len(front)-1] for i := 0; i < len(lists); i++ { if i == j { continue } curr := lists[i] head := curr[0] if tail.Equals(head) { n := make(GeoJSONLineCoordinatesZ, len(front)+len(curr)-1) copy(n, front[:len(front)-1]) n[len(front)-1] = head.mid(tail) copy(n[len(front):], curr[1:]) lists[j] = n if i < len(lists)-1 { copy(lists[i:], lists[i+1:]) } lists[len(lists)-1] = nil lists = lists[:len(lists)-1] continue next } } j++ } return lists } */