Mercurial > gemma
view pkg/imports/sr.go @ 964:1e2dce348cfb
Serialize boundary polygon of sounding result as WKB.
author | Sascha L. Teichmann <sascha.teichmann@intevation.de> |
---|---|
date | Wed, 17 Oct 2018 15:29:58 +0200 |
parents | a4c92e0ef2e1 |
children | f9fb6c399f3f |
line wrap: on
line source
package imports import ( "archive/zip" "bufio" "bytes" "context" "database/sql" "encoding/binary" "encoding/json" "errors" "fmt" "io" "log" "math" "os" "path" "strconv" "strings" "time" shp "github.com/jonas-p/go-shp" "gemma.intevation.de/gemma/pkg/octree" ) type SoundingResult struct { who string zip string } const SoundingResultDateFormat = "2006-01-02" type ( SoundingResultDate struct{ time.Time } SoundingResultMeta struct { Date SoundingResultDate `json:"date"` Bottleneck string `json:"bottleneck"` EPSG uint `json:"epsg"` DepthReference string `json:"depth-reference"` } Point struct { X float64 Y float64 } LineString []Point Polygon []LineString ) const ( wkbNDR byte = 1 wkbPolygon uint32 = 3 ) const ( insertPointsSQL = ` INSERT INTO waterway.sounding_results ( bottleneck_id, date_info, depth_reference, point_cloud, area ) VALUES ( (SELECT bottleneck_id from waterway.bottlenecks where objnam = $1), $2::date, $3, ST_Transform(ST_GeomFromWKB($4, $6::integer), 4326)::geography, (SELECT CASE $5 IS NULL THEN ST_Transform( ST_ConcaveHull( ST_Force2D(ST_GeomFromWKB($4, $6::integer)), 0.7), 4326)::geography ELSE ST_Transform(ST_GeomFromWKB($5, $6::integer), 4326)::geography END) ) RETURNING id` ) func (srd *SoundingResultDate) UnmarshalJSON(data []byte) error { var s string if err := json.Unmarshal(data, &s); err != nil { return err } d, err := time.Parse(SoundingResultDateFormat, s) if err == nil { *srd = SoundingResultDate{d} } return err } func (sr *SoundingResult) Who() string { return sr.who } func (sr *SoundingResult) CleanUp() error { return os.RemoveAll(sr.zip) } func find(needle string, haystack []*zip.File) *zip.File { needle = strings.ToLower(needle) for _, straw := range haystack { if strings.HasSuffix(strings.ToLower(straw.Name), needle) { return straw } } return nil } func loadMeta(f *zip.File) (*SoundingResultMeta, error) { r, err := f.Open() if err != nil { return nil, err } defer r.Close() var m SoundingResultMeta err = json.NewDecoder(r).Decode(&m) return &m, err } func (m *SoundingResultMeta) validate(conn *sql.Conn) error { var b bool err := conn.QueryRowContext(context.Background(), `SELECT true FROM internal.depth_references WHERE depth_reference = $1`, m.DepthReference).Scan(&b) switch { case err == sql.ErrNoRows: return fmt.Errorf("Unknown depth reference '%s'\n", m.DepthReference) case err != nil: return err case !b: return errors.New("Unexpected depth reference") } err = conn.QueryRowContext(context.Background(), `SELECT true FROM waterway.bottlenecks WHERE objnam = $1`, m.Bottleneck).Scan(&b) switch { case err == sql.ErrNoRows: return fmt.Errorf("Unknown bottleneck '%s'\n", m.Bottleneck) case err != nil: return err case !b: return errors.New("Unexpected bottleneck") } return nil } func loadXYZReader(r io.Reader) (octree.MultiPointZ, error) { mpz := make(octree.MultiPointZ, 0, 250000) s := bufio.NewScanner(r) for line := 1; s.Scan(); line++ { text := s.Text() var p octree.Vertex // fmt.Sscanf(text, "%f,%f,%f") is 4 times slower. idx := strings.IndexByte(text, ',') if idx == -1 { log.Printf("format error in line %d\n", line) continue } var err error if p.X, err = strconv.ParseFloat(text[:idx], 64); err != nil { log.Printf("format error in line %d: %v\n", line, err) continue } text = text[idx+1:] if idx = strings.IndexByte(text, ','); idx == -1 { log.Printf("format error in line %d\n", line) continue } if p.Y, err = strconv.ParseFloat(text[:idx], 64); err != nil { log.Printf("format error in line %d: %v\n", line, err) continue } text = text[idx+1:] if p.Z, err = strconv.ParseFloat(text, 64); err != nil { log.Printf("format error in line %d: %v\n", line, err) continue } mpz = append(mpz, p) } if err := s.Err(); err != nil { return nil, err } return mpz, nil } func loadXYZ(f *zip.File) (octree.MultiPointZ, error) { r, err := f.Open() if err != nil { return nil, err } defer r.Close() return loadXYZReader(r) } func toPolygon(numParts int32, parts []int32, points []shp.Point) Polygon { out := make(Polygon, numParts) pos := 0 for i := range out { ps := parts[i] line := make(LineString, ps) for j := int32(0); j < ps; j, pos = j+1, pos+1 { p := &points[pos] line[j] = Point{p.X, p.Y} } out[i] = line } return out } func loadBoundary(z *zip.ReadCloser) (Polygon, error) { shpF := find(".shp", z.File) if shpF == nil { return nil, nil } prefix := strings.TrimSuffix(shpF.Name, path.Ext(shpF.Name)) dbfF := find(prefix+".dbf", z.File) if dbfF == nil { return nil, fmt.Errorf("No DBF file found for %s", shpF.Name) } shpR, err := shpF.Open() if err != nil { return nil, err } dbfR, err := dbfF.Open() if err != nil { shpR.Close() return nil, err } sr := shp.SequentialReaderFromExt(shpR, dbfR) defer sr.Close() if !sr.Next() { return nil, sr.Err() } _, s := sr.Shape() if s == nil { return nil, sr.Err() } switch p := s.(type) { case *shp.Polygon: return toPolygon(p.NumParts, p.Parts, p.Points), nil case *shp.PolygonZ: return toPolygon(p.NumParts, p.Parts, p.Points), nil case *shp.PolygonM: return toPolygon(p.NumParts, p.Parts, p.Points), nil } return nil, fmt.Errorf("Unsupported shape type %T", s) } func (sr *SoundingResult) Do(conn *sql.Conn) error { z, err := zip.OpenReader(sr.zip) if err != nil { return err } defer z.Close() mf := find("meta.json", z.File) if mf == nil { return errors.New("Cannot find 'meta.json'") } m, err := loadMeta(mf) if err != nil { return err } if err := m.validate(conn); err != nil { return err } xyzf := find(".xyz", z.File) if xyzf == nil { return errors.New("Cannot find any *.xyz file") } xyz, err := loadXYZ(xyzf) if err != nil { return err } if len(xyz) == 0 { return errors.New("XYZ does not contain any vertices.") } // Is there a boundary shapefile in the ZIP archive? polygon, err := loadBoundary(z) if err != nil { return err } tx, err := conn.BeginTx(context.Background(), nil) if err != nil { return err } defer tx.Rollback() var id int64 err = tx.QueryRow(insertPointsSQL, m.Bottleneck, m.Date.Time, m.DepthReference, xyz.AsWKB(), polygon.AsWBK(), m.EPSG).Scan(&id) if err != nil { return err } // TODO: Build octree // TODO: Generate iso-lines return tx.Commit() } func (p Polygon) AsWBK() []byte { if p == nil { return nil } // pre-calculate size to avoid reallocations. size := 1 + 4 + 4 for _, ring := range p { size += 4 + len(ring)*2*8 } buf := bytes.NewBuffer(make([]byte, 0, size)) binary.Write(buf, binary.LittleEndian, wkbNDR) binary.Write(buf, binary.LittleEndian, wkbPolygon) binary.Write(buf, binary.LittleEndian, uint32(len(p))) for _, ring := range p { binary.Write(buf, binary.LittleEndian, uint32(len(ring))) for _, v := range ring { binary.Write(buf, binary.LittleEndian, math.Float64bits(v.X)) binary.Write(buf, binary.LittleEndian, math.Float64bits(v.Y)) } } return buf.Bytes() }