From 731cf88dcfdef6b0cf5060501cd753d0828df78f Mon Sep 17 00:00:00 2001 From: Shane Date: Tue, 29 Dec 2020 16:34:08 -0500 Subject: [PATCH 01/14] Extract storage area polygon --- tools/geom.go | 48 +++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 47 insertions(+), 1 deletion(-) diff --git a/tools/geom.go b/tools/geom.go index 17be292..265eca3 100644 --- a/tools/geom.go +++ b/tools/geom.go @@ -100,7 +100,7 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error } // sourceCRS -var sourceCRS string = `PROJCS["NAD_1983_StatePlane_Texas_Central_FIPS_4203_Feet",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["False_Easting",2296583.333333333],PARAMETER["False_Northing",9842500.0],PARAMETER["Central_Meridian",-100.3333333333333],PARAMETER["Standard_Parallel_1",30.11666666666667],PARAMETER["Standard_Parallel_2",31.88333333333333],PARAMETER["Latitude_Of_Origin",29.66666666666667],UNIT["Foot_US",0.3048006096012192]]` +var sourceCRS string = `PROJCS["NAD_1983_StatePlane_Maryland_FIPS_1900_Feet",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["False_Easting",1312333.333333333],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-77.0],PARAMETER["Standard_Parallel_1",38.3],PARAMETER["Standard_Parallel_2",39.45],PARAMETER["Latitude_Of_Origin",37.66666666666666],UNIT["Foot_US",0.3048006096012192]]` // DestinationCRS ... var DestinationCRS int = 4326 @@ -283,6 +283,17 @@ func flipXYLineString25D(xyzLineString gdal.Geometry) gdal.Geometry { return yxzLineString } +func flipXYLinearRing(xyLinearRing gdal.Geometry) gdal.Geometry { + yxLinearRing := gdal.Create(gdal.GT_LinearRing) + nPoints := xyLinearRing.PointCount() + for i := 0; i < nPoints; i++ { + x, y, _ := xyLinearRing.Point(i) + yxLinearRing.AddPoint2D(y, x) + } + xyLinearRing.Destroy() + return yxLinearRing +} + func flipXYPoint(xyPoint gdal.Geometry) gdal.Geometry { yxPoint := gdal.Create(gdal.GT_Point) nPoints := xyPoint.PointCount() @@ -405,6 +416,34 @@ func getBanks(line string, transform gdal.CoordinateTransform, xsLayer vectorLay return layers, nil } +func getStorageArea(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vectorLayer, error) { + layer := vectorLayer{FeatureName: strings.TrimSpace(strings.Split(rightofEquals(sc.Text()), ",")[0])} + + xyPairs, err := getDataPairsfromTextBlock("Storage Area Surface Line=", sc, 32, 16) + if err != nil { + return layer, err + } + + xyLinearRing := gdal.Create(gdal.GT_LinearRing) + for _, pair := range xyPairs { + xyLinearRing.AddPoint2D(pair[0], pair[1]) + } + + xyLinearRing.Transform(transform) + // This is a temporary fix since the x and y values need to be flipped: + yxLinearRing := flipXYLinearRing(xyLinearRing) + + yxPolygon := gdal.Create(gdal.GT_Polygon) + yxPolygon.AddGeometry(yxLinearRing) + yxMultiPolygon := yxPolygon.ForceToMultiPolygon() + wkb, err := yxMultiPolygon.ToWKB() + if err != nil { + return layer, err + } + layer.Geometry = wkb + return layer, err +} + // GetGeospatialData ... func GetGeospatialData(gd *GeoData, fs filestore.FileStore, geomFilePath string) error { geomFileName := filepath.Base(geomFilePath) @@ -435,6 +474,13 @@ func GetGeospatialData(gd *GeoData, fs filestore.FileStore, geomFilePath string) f.Rivers = append(f.Rivers, riverLayer) riverReachName = riverLayer.FeatureName + case strings.HasPrefix(line, "Storage Area="): + storageAreaLayer, err := getStorageArea(sc, transform) + if err != nil { + return err + } + f.StorageAreas = append(f.StorageAreas, storageAreaLayer) + case strings.HasPrefix(line, "Type RM Length L Ch R = 1"): xsLayer, bankLayers, err := getXSBanks(sc, transform, riverReachName) if err != nil { From 5e79d90fbb738a5520757032e1f0d7161d8f2061 Mon Sep 17 00:00:00 2001 From: Shane Date: Wed, 30 Dec 2020 10:23:44 -0500 Subject: [PATCH 02/14] Count the number of hydraulic structures --- handlers/index.go | 10 +-- tools/geom.go | 151 ++++++++++++++++++++++++++++++---------------- 2 files changed, 104 insertions(+), 57 deletions(-) diff --git a/handlers/index.go b/handlers/index.go index e994be6..79e3f2e 100644 --- a/handlers/index.go +++ b/handlers/index.go @@ -28,11 +28,11 @@ func Index(fs *filestore.FileStore) echo.HandlerFunc { if err != nil { return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) } - mod, err := rm.Index() - if err != nil { - return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) - } + // mod, err := rm.Index() + // if err != nil { + // return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) + // } - return c.JSON(http.StatusOK, mod) + return c.JSON(http.StatusOK, rm) } } diff --git a/tools/geom.go b/tools/geom.go index 265eca3..8c40dd8 100644 --- a/tools/geom.go +++ b/tools/geom.go @@ -18,10 +18,83 @@ import ( // GeomFileContents keywords and data container for ras flow file search type GeomFileContents struct { Path string - FileExt string //`json:"File Extension"` - GeomTitle string //`json:"Geom Title"` - ProgramVersion string //`json:"Program Version"` - Description string //`json:"Description"` + FileExt string //`json:"File Extension"` + GeomTitle string //`json:"Geom Title"` + ProgramVersion string //`json:"Program Version"` + Description string //`json:"Description"` + Structures []hydraulicStructures //`json:"Hydraulic Structures"` +} + +type hydraulicStructures struct { + River string //`json:"River Name"` + Reach string //`json:"Reach Name"` + NumXS int //`json:"N CrossSections"` + NumCulverts int //`json:"N Culverts"` + BridgeData bridgeData //`json:"Bridges"` + NumInlines int //`json:"N Inlines"` +} + +type bridgeData struct { + NumBridges int //`json:"Num Bridges"` + Bridges []bridges //`json:"Num Bridges"` +} + +type bridges struct { + Name string +} + +func rightofEquals(line string) string { + return strings.TrimSpace(strings.Split(line, "=")[1]) +} + +func getDescription(sc *bufio.Scanner) (string, error) { + description := "" + for sc.Scan() { + line := sc.Text() + endDescription, err := regexp.MatchString("END GEOM DESCRIPTION", line) + if err != nil { + return description, err + } + if endDescription { + return description, nil + } + if line != "" { + description += line + "\n" + } + } + return description, nil +} + +func getHydraulicStructureData(sc *bufio.Scanner) (hydraulicStructures, error) { + riverReach := strings.Split(rightofEquals(sc.Text()), ",") + structures := hydraulicStructures{River: strings.TrimSpace(riverReach[0]), Reach: strings.TrimSpace(riverReach[1])} + bData := bridgeData{} + for sc.Scan() { + line := sc.Text() + if strings.HasPrefix(line, "Type RM Length L Ch R =") { + data := strings.Split(rightofEquals(line), ",") + structureType, err := strconv.Atoi(strings.TrimSpace(data[0])) + if err != nil { + return structures, err + } + switch structureType { + case 1: + structures.NumXS++ + + case 2: + structures.NumCulverts++ + + case 3: + bData.NumBridges++ + + case 5: + structures.NumInlines++ + + } + } + } + structures.BridgeData = bData + return structures, nil } // getGeomData Reads a geometry file. returns none to allow concurrency @@ -42,59 +115,37 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error var line string header := true for sc.Scan() { - line = sc.Text() + switch { + case strings.HasPrefix(line, "Geom Title="): + meta.GeomTitle = rightofEquals(line) - match, err := regexp.MatchString("=", line) - - if err != nil { - errChan <- err - return - } - - beginDescription, err := regexp.MatchString("BEGIN GEOM DESCRIPTION", line) - - if err != nil { - errChan <- err - return - } - - if match { - data := strings.Split(line, "=") - - switch data[0] { - - case "Geom Title": - meta.GeomTitle = data[1] - - case "Program Version": - meta.ProgramVersion = data[1] - - case "River Reach", "Storage Area": - header = false - - } - - } else if header && beginDescription { - - for sc.Scan() { - line = sc.Text() - endDescription, _ := regexp.MatchString("END GEOM DESCRIPTION", line) - - if endDescription { - break + case strings.HasPrefix(line, "Program Version="): + meta.ProgramVersion = rightofEquals(line) - } else { - if line != "" { - meta.Description += line + "\n" - } + case strings.HasPrefix(line, "BEGIN GEOM DESCRIPTION:"): + if header { + description, err := getDescription(sc) + if err != nil { + errChan <- err + return } + meta.Description += description + } + case strings.HasPrefix(line, "River Reach="): + structures, err := getHydraulicStructureData(sc) + if err != nil { + errChan <- err + return } + meta.Structures = append(meta.Structures, structures) + header = false + case strings.HasPrefix(line, "Storage Area="): + header = false } } - rm.Metadata.GeomFiles = append(rm.Metadata.GeomFiles, meta) return } @@ -133,10 +184,6 @@ type xyzPoint struct { z float64 } -func rightofEquals(line string) string { - return strings.TrimSpace(strings.Split(line, "=")[1]) -} - func dataPairsfromTextBlock(sc *bufio.Scanner, nPairs int, colWidth int, valueWidth int) ([][2]float64, error) { var stride int = valueWidth * 2 pairs := [][2]float64{} From d3c57147627fa5e71456511124442d161e2cf44a Mon Sep 17 00:00:00 2001 From: Shane Date: Wed, 30 Dec 2020 20:40:09 -0500 Subject: [PATCH 03/14] Extract bridge metadata from geometry file --- tools/geom.go | 231 +++++++++++++++++++++++++++++++++++++++++++++----- 1 file changed, 212 insertions(+), 19 deletions(-) diff --git a/tools/geom.go b/tools/geom.go index 8c40dd8..33a8e91 100644 --- a/tools/geom.go +++ b/tools/geom.go @@ -18,40 +18,84 @@ import ( // GeomFileContents keywords and data container for ras flow file search type GeomFileContents struct { Path string - FileExt string //`json:"File Extension"` - GeomTitle string //`json:"Geom Title"` - ProgramVersion string //`json:"Program Version"` - Description string //`json:"Description"` - Structures []hydraulicStructures //`json:"Hydraulic Structures"` + FileExt string `json:"File Extension"` + GeomTitle string `json:"Geom Title"` + ProgramVersion string `json:"Program Version"` + Description string `json:"Description"` + Structures []hydraulicStructures `json:"Hydraulic Structures"` } type hydraulicStructures struct { - River string //`json:"River Name"` - Reach string //`json:"Reach Name"` - NumXS int //`json:"N CrossSections"` - NumCulverts int //`json:"N Culverts"` - BridgeData bridgeData //`json:"Bridges"` - NumInlines int //`json:"N Inlines"` + River string `json:"River Name"` + Reach string `json:"Reach Name"` + NumXS int `json:"Num CrossSections"` + NumCulverts int `json:"Num Culverts"` + BridgeData bridgeData `json:"Bridges"` + NumInlines int `json:"Num Inlines"` } type bridgeData struct { - NumBridges int //`json:"Num Bridges"` - Bridges []bridges //`json:"Num Bridges"` + NumBridges int `json:"Num Bridges"` + Bridges []bridges `json:"Bridges"` } type bridges struct { - Name string + Name string + Station int + Description string + DeckWidth float64 `json:"Deck Width"` + UpHighChord chordPairs `json:"Upstream High Chord"` + UpLowChord chordPairs `json:"Upstream Low Chord"` + DownHighChord chordPairs `json:"Downstream High Chord"` + DownLowChord chordPairs `json:"Downstream Max Chord"` + NumPiers int `json:"Num Piers"` +} + +type chordPairs struct { + Max float64 + Min float64 +} + +func minValue(values []float64) (float64, error) { + if len(values) == 0 { + return 0.0, errors.New("Cannot detect a minimum value in an empty slice") + } + + min := values[0] + for _, v := range values { + if v < min { + min = v + } + } + + return min, nil +} + +func maxValue(values []float64) (float64, error) { + if len(values) == 0 { + return 0.0, errors.New("Cannot detect a maximum value in an empty slice") + } + + max := values[0] + for _, v := range values { + if v > max { + max = v + } + } + + return max, nil } func rightofEquals(line string) string { return strings.TrimSpace(strings.Split(line, "=")[1]) } -func getDescription(sc *bufio.Scanner) (string, error) { +func getDescription(sc *bufio.Scanner, endLine string) (string, error) { description := "" + nLines := 0 for sc.Scan() { line := sc.Text() - endDescription, err := regexp.MatchString("END GEOM DESCRIPTION", line) + endDescription, err := regexp.MatchString(endLine, line) if err != nil { return description, err } @@ -59,11 +103,155 @@ func getDescription(sc *bufio.Scanner) (string, error) { return description, nil } if line != "" { - description += line + "\n" + if nLines > 0 { + description += "\n" + } + description += line + nLines++ } } return description, nil } +func numberofLines(nValues int, colWidth int, valueWidth int) int { + nLines := math.Ceil(float64(nValues) / (float64(colWidth) / float64(valueWidth))) + return int(nLines) +} + +func datafromTextBlock(sc *bufio.Scanner, colWidth int, valueWidth int, nLines int, nSkipLines int) ([]float64, error) { + values := []float64{} + nSkipped := 0 + nProcessed := 0 +out: + for sc.Scan() { + if nSkipped < nSkipLines { + nSkipped++ + continue + } + nProcessed++ + line := sc.Text() + for s := 0; s < colWidth; { + if len(line) > s { + sVal := strings.TrimSpace(line[s : s+valueWidth]) + if sVal != "" { + val, err := strconv.ParseFloat(sVal, 64) + if err != nil { + return values, err + } + values = append(values, val) + } + s += valueWidth + } else { + if nLines == nProcessed { + break out + } + break + } + + } + if nLines == nProcessed { + break out + } + } + return values, nil +} + +func getHighLowChord(sc *bufio.Scanner, nElevsText string, colWidth int, valueWidth int) ([2]chordPairs, error) { + highLowPairs := [2]chordPairs{} + + nElevs, err := strconv.Atoi(strings.TrimSpace(nElevsText)) + if err != nil { + return highLowPairs, err + } + + nLines := numberofLines(nElevs, colWidth, valueWidth) + + elevHighChord, err := datafromTextBlock(sc, colWidth, valueWidth, nLines, nLines) + if err != nil { + return highLowPairs, err + } + + maxHighCord, err := maxValue(elevHighChord) + if err != nil { + return highLowPairs, err + } + + minHighCord, err := minValue(elevHighChord) + if err != nil { + return highLowPairs, err + } + highLowPairs[0] = chordPairs{Max: maxHighCord, Min: minHighCord} + + elevLowChord, err := datafromTextBlock(sc, 80, 8, nLines, 0) + if err != nil { + return highLowPairs, err + } + + maxLowCord, err := maxValue(elevLowChord) + if err != nil { + return highLowPairs, err + } + + minLowCord, err := minValue(elevLowChord) + if err != nil { + return highLowPairs, err + } + highLowPairs[1] = chordPairs{Max: maxLowCord, Min: minLowCord} + + return highLowPairs, nil +} + +func getBridgeData(sc *bufio.Scanner, lineData []string) (bridges, error) { + bridge := bridges{} + station, err := strconv.Atoi(strings.TrimSpace(lineData[1])) + if err != nil { + return bridge, err + } + bridge.Station = station + + for sc.Scan() { + line := sc.Text() + switch { + case strings.HasPrefix(line, "BEGIN DESCRIPTION"): + description, err := getDescription(sc, "END DESCRIPTION:") + if err != nil { + return bridge, err + } + bridge.Description += description + + case strings.HasPrefix(line, "Node Name="): + bridge.Name = rightofEquals(line) + + case strings.HasPrefix(line, "Deck Dist"): + sc.Scan() + nextLineData := strings.Split(sc.Text(), ",") + deckWidth, err := strconv.ParseFloat(strings.TrimSpace(nextLineData[0]), 64) + if err != nil { + return bridge, err + } + bridge.DeckWidth = deckWidth + upHighLowPair, err := getHighLowChord(sc, nextLineData[4], 80, 8) + if err != nil { + return bridge, err + } + bridge.UpHighChord = upHighLowPair[0] + bridge.UpLowChord = upHighLowPair[1] + + downHighLowPair, err := getHighLowChord(sc, nextLineData[5], 80, 8) + if err != nil { + return bridge, err + } + bridge.DownHighChord = downHighLowPair[0] + bridge.DownLowChord = downHighLowPair[1] + + case strings.HasPrefix(line, "Pier Skew"): + bridge.NumPiers++ + + case strings.HasPrefix(line, "BR Coef"): + return bridge, err + } + } + return bridge, nil +} func getHydraulicStructureData(sc *bufio.Scanner) (hydraulicStructures, error) { riverReach := strings.Split(rightofEquals(sc.Text()), ",") @@ -85,6 +273,11 @@ func getHydraulicStructureData(sc *bufio.Scanner) (hydraulicStructures, error) { structures.NumCulverts++ case 3: + bridge, err := getBridgeData(sc, data) + if err != nil { + return structures, err + } + bData.Bridges = append(bData.Bridges, bridge) bData.NumBridges++ case 5: @@ -125,7 +318,7 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error case strings.HasPrefix(line, "BEGIN GEOM DESCRIPTION:"): if header { - description, err := getDescription(sc) + description, err := getDescription(sc, "END GEOM DESCRIPTION:") if err != nil { errChan <- err return @@ -201,7 +394,7 @@ out: return pairs, err } pairs = append(pairs, [2]float64{val1, val2}) - if int(len(pairs)) == nPairs { + if len(pairs) == nPairs { break out } } else { From 3b15e07d9bf53e40c24d782f581e79f4704ec3ff Mon Sep 17 00:00:00 2001 From: Shane Date: Thu, 31 Dec 2020 10:06:26 -0500 Subject: [PATCH 04/14] Extract bridges when there are multiple river reaches --- tools/geom.go | 146 +++++++++++++++++++++++++++++--------------------- 1 file changed, 86 insertions(+), 60 deletions(-) diff --git a/tools/geom.go b/tools/geom.go index 33a8e91..e9a61bf 100644 --- a/tools/geom.go +++ b/tools/geom.go @@ -41,7 +41,7 @@ type bridgeData struct { type bridges struct { Name string - Station int + Station float64 Description string DeckWidth float64 `json:"Deck Width"` UpHighChord chordPairs `json:"Upstream High Chord"` @@ -56,51 +56,52 @@ type chordPairs struct { Min float64 } -func minValue(values []float64) (float64, error) { +func maxValue(values []float64) (float64, error) { if len(values) == 0 { - return 0.0, errors.New("Cannot detect a minimum value in an empty slice") + return 0.0, errors.New("Cannot detect a maximum value in an empty slice") } - min := values[0] + max := values[0] for _, v := range values { - if v < min { - min = v + if v > max { + max = v } } - return min, nil + return max, nil } -func maxValue(values []float64) (float64, error) { +func minValue(values []float64) (float64, error) { if len(values) == 0 { - return 0.0, errors.New("Cannot detect a maximum value in an empty slice") + return 0.0, errors.New("Cannot detect a minimum value in an empty slice") } - max := values[0] + min := values[0] for _, v := range values { - if v > max { - max = v + if v < min { + min = v } } - return max, nil + return min, nil } func rightofEquals(line string) string { return strings.TrimSpace(strings.Split(line, "=")[1]) } -func getDescription(sc *bufio.Scanner, endLine string) (string, error) { +func getDescription(sc *bufio.Scanner, idx int, endLine string) (string, int, error) { description := "" nLines := 0 for sc.Scan() { line := sc.Text() + idx++ endDescription, err := regexp.MatchString(endLine, line) if err != nil { - return description, err + return description, idx, err } if endDescription { - return description, nil + return description, idx, nil } if line != "" { if nLines > 0 { @@ -110,25 +111,25 @@ func getDescription(sc *bufio.Scanner, endLine string) (string, error) { nLines++ } } - return description, nil + return description, idx, nil } func numberofLines(nValues int, colWidth int, valueWidth int) int { nLines := math.Ceil(float64(nValues) / (float64(colWidth) / float64(valueWidth))) return int(nLines) } -func datafromTextBlock(sc *bufio.Scanner, colWidth int, valueWidth int, nLines int, nSkipLines int) ([]float64, error) { +func datafromTextBlock(hsSc *bufio.Scanner, colWidth int, valueWidth int, nLines int, nSkipLines int) ([]float64, error) { values := []float64{} nSkipped := 0 nProcessed := 0 out: - for sc.Scan() { + for hsSc.Scan() { if nSkipped < nSkipLines { nSkipped++ continue } nProcessed++ - line := sc.Text() + line := hsSc.Text() for s := 0; s < colWidth; { if len(line) > s { sVal := strings.TrimSpace(line[s : s+valueWidth]) @@ -155,7 +156,7 @@ out: return values, nil } -func getHighLowChord(sc *bufio.Scanner, nElevsText string, colWidth int, valueWidth int) ([2]chordPairs, error) { +func getHighLowChord(hsSc *bufio.Scanner, nElevsText string, colWidth int, valueWidth int) ([2]chordPairs, error) { highLowPairs := [2]chordPairs{} nElevs, err := strconv.Atoi(strings.TrimSpace(nElevsText)) @@ -165,7 +166,7 @@ func getHighLowChord(sc *bufio.Scanner, nElevsText string, colWidth int, valueWi nLines := numberofLines(nElevs, colWidth, valueWidth) - elevHighChord, err := datafromTextBlock(sc, colWidth, valueWidth, nLines, nLines) + elevHighChord, err := datafromTextBlock(hsSc, colWidth, valueWidth, nLines, nLines) if err != nil { return highLowPairs, err } @@ -181,7 +182,7 @@ func getHighLowChord(sc *bufio.Scanner, nElevsText string, colWidth int, valueWi } highLowPairs[0] = chordPairs{Max: maxHighCord, Min: minHighCord} - elevLowChord, err := datafromTextBlock(sc, 80, 8, nLines, 0) + elevLowChord, err := datafromTextBlock(hsSc, 80, 8, nLines, 0) if err != nil { return highLowPairs, err } @@ -196,23 +197,23 @@ func getHighLowChord(sc *bufio.Scanner, nElevsText string, colWidth int, valueWi return highLowPairs, err } highLowPairs[1] = chordPairs{Max: maxLowCord, Min: minLowCord} - return highLowPairs, nil } -func getBridgeData(sc *bufio.Scanner, lineData []string) (bridges, error) { +func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { bridge := bridges{} - station, err := strconv.Atoi(strings.TrimSpace(lineData[1])) + + station, err := strconv.ParseFloat(strings.TrimSpace(lineData[1]), 64) if err != nil { return bridge, err } bridge.Station = station - for sc.Scan() { - line := sc.Text() + for hsSc.Scan() { + line := hsSc.Text() switch { case strings.HasPrefix(line, "BEGIN DESCRIPTION"): - description, err := getDescription(sc, "END DESCRIPTION:") + description, _, err := getDescription(hsSc, 0, "END DESCRIPTION:") if err != nil { return bridge, err } @@ -222,21 +223,21 @@ func getBridgeData(sc *bufio.Scanner, lineData []string) (bridges, error) { bridge.Name = rightofEquals(line) case strings.HasPrefix(line, "Deck Dist"): - sc.Scan() - nextLineData := strings.Split(sc.Text(), ",") + hsSc.Scan() + nextLineData := strings.Split(hsSc.Text(), ",") deckWidth, err := strconv.ParseFloat(strings.TrimSpace(nextLineData[0]), 64) if err != nil { return bridge, err } bridge.DeckWidth = deckWidth - upHighLowPair, err := getHighLowChord(sc, nextLineData[4], 80, 8) + upHighLowPair, err := getHighLowChord(hsSc, nextLineData[4], 80, 8) if err != nil { return bridge, err } bridge.UpHighChord = upHighLowPair[0] bridge.UpLowChord = upHighLowPair[1] - downHighLowPair, err := getHighLowChord(sc, nextLineData[5], 80, 8) + downHighLowPair, err := getHighLowChord(hsSc, nextLineData[5], 80, 8) if err != nil { return bridge, err } @@ -253,40 +254,61 @@ func getBridgeData(sc *bufio.Scanner, lineData []string) (bridges, error) { return bridge, nil } -func getHydraulicStructureData(sc *bufio.Scanner) (hydraulicStructures, error) { - riverReach := strings.Split(rightofEquals(sc.Text()), ",") - structures := hydraulicStructures{River: strings.TrimSpace(riverReach[0]), Reach: strings.TrimSpace(riverReach[1])} +func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStructures, error) { + structures := hydraulicStructures{} bData := bridgeData{} - for sc.Scan() { - line := sc.Text() - if strings.HasPrefix(line, "Type RM Length L Ch R =") { - data := strings.Split(rightofEquals(line), ",") - structureType, err := strconv.Atoi(strings.TrimSpace(data[0])) - if err != nil { - return structures, err - } - switch structureType { - case 1: - structures.NumXS++ - - case 2: - structures.NumCulverts++ - case 3: - bridge, err := getBridgeData(sc, data) + newf, err := rm.FileStore.GetObject(fn) + if err != nil { + return structures, nil + } + defer newf.Close() + + hsSc := bufio.NewScanner(newf) + + i := 0 + for hsSc.Scan() { + if i == idx { + riverReach := strings.Split(rightofEquals(hsSc.Text()), ",") + structures.River = strings.TrimSpace(riverReach[0]) + structures.Reach = strings.TrimSpace(riverReach[1]) + } else if i > idx { + line := hsSc.Text() + if strings.HasPrefix(line, "River Reach=") { + structures.BridgeData = bData + return structures, nil + } + if strings.HasPrefix(line, "Type RM Length L Ch R =") { + data := strings.Split(rightofEquals(line), ",") + structureType, err := strconv.Atoi(strings.TrimSpace(data[0])) if err != nil { return structures, err } - bData.Bridges = append(bData.Bridges, bridge) - bData.NumBridges++ + switch structureType { + case 1: + structures.NumXS++ + + case 2: + structures.NumCulverts++ + + case 3: + bridge, err := getBridgeData(hsSc, data) + if err != nil { + return structures, err + } + bData.Bridges = append(bData.Bridges, bridge) + bData.NumBridges++ - case 5: - structures.NumInlines++ + case 5: + structures.NumInlines++ + } } } + i++ } structures.BridgeData = bData + return structures, nil } @@ -305,10 +327,13 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error defer f.Close() sc := bufio.NewScanner(f) - var line string + + var description string + header := true + idx := 0 for sc.Scan() { - line = sc.Text() + line := sc.Text() switch { case strings.HasPrefix(line, "Geom Title="): meta.GeomTitle = rightofEquals(line) @@ -318,7 +343,7 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error case strings.HasPrefix(line, "BEGIN GEOM DESCRIPTION:"): if header { - description, err := getDescription(sc, "END GEOM DESCRIPTION:") + description, idx, err = getDescription(sc, idx, "END GEOM DESCRIPTION:") if err != nil { errChan <- err return @@ -327,7 +352,7 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error } case strings.HasPrefix(line, "River Reach="): - structures, err := getHydraulicStructureData(sc) + structures, err := getHydraulicStructureData(rm, fn, idx) if err != nil { errChan <- err return @@ -338,6 +363,7 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error case strings.HasPrefix(line, "Storage Area="): header = false } + idx++ } rm.Metadata.GeomFiles = append(rm.Metadata.GeomFiles, meta) return From a106aa95c87cb4c10d1705e22cbabc830d32f0e3 Mon Sep 17 00:00:00 2001 From: Shane Date: Thu, 31 Dec 2020 11:35:31 -0500 Subject: [PATCH 05/14] Read projection file and pass source crs to geospatial extraction function --- tools/geom.go | 35 ++++++++++++++++++++++++++--------- tools/model.go | 26 +++++++++++++++++++++----- tools/project.go | 1 + 3 files changed, 48 insertions(+), 14 deletions(-) diff --git a/tools/geom.go b/tools/geom.go index e9a61bf..401d8fa 100644 --- a/tools/geom.go +++ b/tools/geom.go @@ -369,11 +369,28 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error return } -// sourceCRS -var sourceCRS string = `PROJCS["NAD_1983_StatePlane_Maryland_FIPS_1900_Feet",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["False_Easting",1312333.333333333],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-77.0],PARAMETER["Standard_Parallel_1",38.3],PARAMETER["Standard_Parallel_2",39.45],PARAMETER["Latitude_Of_Origin",37.66666666666666],UNIT["Foot_US",0.3048006096012192]]` +// getProjection Reads a projection file. returns none to allow concurrency +func getProjection(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error) { -// DestinationCRS ... -var DestinationCRS int = 4326 + defer wg.Done() + + f, err := rm.FileStore.GetObject(fn) + if err != nil { + errChan <- err + return + } + defer f.Close() + sc := bufio.NewScanner(f) + sc.Scan() + line := sc.Text() + if strings.HasPrefix(line, "PROJCS[") { + rm.Metadata.Projection = line + return + } + + errChan <- errors.New("Unexpected projection file structure") + return +} // GeoData ... type GeoData struct { @@ -511,16 +528,16 @@ func attributeZ(xyPairs [][2]float64, mzPairs [][2]float64) []xyzPoint { return points } -func getTransform(sourceProjection string, destinationEPSG int) (gdal.CoordinateTransform, error) { +func getTransform(sourceCRS string, destinationCRS int) (gdal.CoordinateTransform, error) { transform := gdal.CoordinateTransform{} - sourceSpRef := gdal.CreateSpatialReference(sourceProjection) + sourceSpRef := gdal.CreateSpatialReference(sourceCRS) sourceSpRef.MorphFromESRI() if err := sourceSpRef.Validate(); err != nil { return transform, errors.New("Unable to extract geospatial data. Invalid source Projection") } destinationSpRef := gdal.CreateSpatialReference("") - if err := destinationSpRef.FromEPSG(destinationEPSG); err != nil { + if err := destinationSpRef.FromEPSG(destinationCRS); err != nil { return transform, err } transform = gdal.CreateCoordinateTransform(sourceSpRef, destinationSpRef) @@ -711,7 +728,7 @@ func getStorageArea(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vect } // GetGeospatialData ... -func GetGeospatialData(gd *GeoData, fs filestore.FileStore, geomFilePath string) error { +func GetGeospatialData(gd *GeoData, fs filestore.FileStore, geomFilePath string, sourceCRS string, destinationCRS int) error { geomFileName := filepath.Base(geomFilePath) f := Features{} riverReachName := "" @@ -724,7 +741,7 @@ func GetGeospatialData(gd *GeoData, fs filestore.FileStore, geomFilePath string) sc := bufio.NewScanner(file) - transform, err := getTransform(sourceCRS, DestinationCRS) + transform, err := getTransform(sourceCRS, destinationCRS) if err != nil { return err } diff --git a/tools/model.go b/tools/model.go index dd707c6..6563e62 100644 --- a/tools/model.go +++ b/tools/model.go @@ -21,6 +21,7 @@ type fileExtMatchers struct { SteadyRun *regexp.Regexp UnsteadyRun *regexp.Regexp AllFlowRun *regexp.Regexp + Projection *regexp.Regexp } var rasRE fileExtMatchers = fileExtMatchers{ // Maybe these ones are better? need a regex experts opinion @@ -34,13 +35,15 @@ var rasRE fileExtMatchers = fileExtMatchers{ // Maybe these ones are better? nee SteadyRun: regexp.MustCompile(".r[0-9][0-9]"), // `^\.r(0[1-9]|[1-9][0-9])$` UnsteadyRun: regexp.MustCompile(".x[0-9][0-9]"), // `^\.x(0[1-9]|[1-9][0-9])$` AllFlowRun: regexp.MustCompile(".[rx][0-9][0-9]"), // `^\.[rx](0[1-9]|[1-9][0-9])$` + Projection: regexp.MustCompile(".projection"), } // holder of multiple wait groups to help process files concurrency type rasWaitGroup struct { - Geom sync.WaitGroup - Plan sync.WaitGroup - Flow sync.WaitGroup + Geom sync.WaitGroup + Plan sync.WaitGroup + Flow sync.WaitGroup + Projection sync.WaitGroup } // Model is a general type should contain all necessary data for a model of any type. @@ -188,12 +191,15 @@ func (rm *RasModel) Index() (Model, error) { return mod, nil } +var destinationCRS int = 4326 + // GeospatialData ... func (rm *RasModel) GeospatialData() (GeoData, error) { - gd := GeoData{Features: make(map[string]Features), Georeference: DestinationCRS} + gd := GeoData{Features: make(map[string]Features), Georeference: destinationCRS} + sourceCRS := rm.Metadata.Projection for _, g := range rm.Metadata.GeomFiles { - if err := GetGeospatialData(&gd, rm.FileStore, g.Path); err != nil { + if err := GetGeospatialData(&gd, rm.FileStore, g.Path, sourceCRS, destinationCRS); err != nil { return gd, err } } @@ -238,6 +244,7 @@ func NewRasModel(key string, fs filestore.FileStore) (*RasModel, error) { errChan := make(chan error) var rasWG rasWaitGroup + nProjection := 0 for _, fp := range rm.FileList { ext := filepath.Ext(fp) @@ -256,12 +263,21 @@ func NewRasModel(key string, fs filestore.FileStore) (*RasModel, error) { rasWG.Flow.Add(1) go getFlowData(&rm, fp, &rasWG.Flow, errChan) + case rasRE.Projection.MatchString(ext): + if nProjection > 0 { + return &rm, errors.New("multiple projection files identified") + } + rasWG.Projection.Add(1) + go getProjection(&rm, fp, &rasWG.Projection, errChan) + nProjection++ + } } rasWG.Plan.Wait() rasWG.Geom.Wait() rasWG.Flow.Wait() + rasWG.Projection.Wait() if len(errChan) > 0 { fmt.Printf("Encountered %d errors\n", len(errChan)) diff --git a/tools/project.go b/tools/project.go index b2567cc..d60402c 100644 --- a/tools/project.go +++ b/tools/project.go @@ -17,6 +17,7 @@ type ProjectMetadata struct { PlanFiles []PlanFileContents //`json:"Plan Data"` FlowFiles []FlowFileContents //`json:"Flow Data"` GeomFiles []GeomFileContents //`json:"Geometry Data"` + Projection string //`json:"Projection"` Notes string //`json:"Notes"` } From 347ef58a6bbf33d621d84c3bad3f9da901100ff0 Mon Sep 17 00:00:00 2001 From: Shane Date: Mon, 4 Jan 2021 08:45:24 -0500 Subject: [PATCH 06/14] Move destination CRS to config file, resolves #4 --- config/config.go | 8 +++++--- handlers/geospatialdata.go | 8 ++++---- main.go | 2 +- tools/model.go | 4 +--- 4 files changed, 11 insertions(+), 11 deletions(-) diff --git a/config/config.go b/config/config.go index b61f6a1..9199692 100644 --- a/config/config.go +++ b/config/config.go @@ -8,9 +8,10 @@ import ( ) type APIConfig struct { - Host string - Port int - FileStore *filestore.FileStore + Host string + Port int + FileStore *filestore.FileStore + DestinationCRS int } // Address tells the application where to run the api out of @@ -24,6 +25,7 @@ func Init() *APIConfig { config.Host = "" // 0.0.0.0 config.Port = 5600 config.FileStore = FileStoreInit(false) + config.DestinationCRS = 4326 return config } diff --git a/handlers/geospatialdata.go b/handlers/geospatialdata.go index 8bb581d..f2e5784 100644 --- a/handlers/geospatialdata.go +++ b/handlers/geospatialdata.go @@ -3,9 +3,9 @@ package handlers import ( "net/http" + "github.com/USACE/mcat-ras/config" ras "github.com/USACE/mcat-ras/tools" - "github.com/USACE/filestore" "github.com/labstack/echo/v4" ) @@ -19,17 +19,17 @@ import ( // @Success 200 {object} interface{} // @Failure 500 {object} SimpleResponse // @Router /geospatialdata [get] -func GeospatialData(fs *filestore.FileStore) echo.HandlerFunc { +func GeospatialData(cf *config.APIConfig) echo.HandlerFunc { return func(c echo.Context) error { definitionFile := c.QueryParam("definition_file") - rm, err := ras.NewRasModel(definitionFile, *fs) + rm, err := ras.NewRasModel(definitionFile, *cf.FileStore) if err != nil { return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) } - data, err := rm.GeospatialData() + data, err := rm.GeospatialData(cf.DestinationCRS) if err != nil { return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) } diff --git a/main.go b/main.go index 743bc50..d413225 100644 --- a/main.go +++ b/main.go @@ -39,7 +39,7 @@ func main() { e.GET("/isgeospatial", handlers.IsGeospatial(appConfig.FileStore)) e.GET("/modeltype", handlers.ModelType(appConfig.FileStore)) e.GET("/modelversion", handlers.ModelVersion(appConfig.FileStore)) - e.GET("/geospatialdata", handlers.GeospatialData(appConfig.FileStore)) + e.GET("/geospatialdata", handlers.GeospatialData(appConfig)) e.Logger.Fatal(e.Start(appConfig.Address())) } diff --git a/tools/model.go b/tools/model.go index 6d218ad..63a9954 100644 --- a/tools/model.go +++ b/tools/model.go @@ -191,10 +191,8 @@ func (rm *RasModel) Index() (Model, error) { return mod, nil } -var destinationCRS int = 4326 - // GeospatialData ... -func (rm *RasModel) GeospatialData() (GeoData, error) { +func (rm *RasModel) GeospatialData(destinationCRS int) (GeoData, error) { gd := GeoData{Features: make(map[string]Features), Georeference: destinationCRS} sourceCRS := rm.Metadata.Projection From 6f843f13b2f83f8249696bc32b1f17b2c0a4d2f3 Mon Sep 17 00:00:00 2001 From: Shane Date: Mon, 4 Jan 2021 10:42:07 -0500 Subject: [PATCH 07/14] Break up the geom.go file, resolves #5 --- handlers/geospatialdata.go | 8 +- main.go | 2 +- tools/geom.go | 706 ------------------------------------- tools/geospatial.go | 401 +++++++++++++++++++++ tools/model.go | 35 +- tools/structures.go | 238 +++++++++++++ tools/utils.go | 66 ++++ 7 files changed, 744 insertions(+), 712 deletions(-) create mode 100644 tools/geospatial.go create mode 100644 tools/structures.go create mode 100644 tools/utils.go diff --git a/handlers/geospatialdata.go b/handlers/geospatialdata.go index f2e5784..505e46b 100644 --- a/handlers/geospatialdata.go +++ b/handlers/geospatialdata.go @@ -3,7 +3,7 @@ package handlers import ( "net/http" - "github.com/USACE/mcat-ras/config" + "github.com/USACE/filestore" ras "github.com/USACE/mcat-ras/tools" "github.com/labstack/echo/v4" @@ -19,17 +19,17 @@ import ( // @Success 200 {object} interface{} // @Failure 500 {object} SimpleResponse // @Router /geospatialdata [get] -func GeospatialData(cf *config.APIConfig) echo.HandlerFunc { +func GeospatialData(fs *filestore.FileStore, destinationCRS int) echo.HandlerFunc { return func(c echo.Context) error { definitionFile := c.QueryParam("definition_file") - rm, err := ras.NewRasModel(definitionFile, *cf.FileStore) + rm, err := ras.NewRasModel(definitionFile, *fs) if err != nil { return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) } - data, err := rm.GeospatialData(cf.DestinationCRS) + data, err := rm.GeospatialData(destinationCRS) if err != nil { return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) } diff --git a/main.go b/main.go index d413225..b6b0361 100644 --- a/main.go +++ b/main.go @@ -39,7 +39,7 @@ func main() { e.GET("/isgeospatial", handlers.IsGeospatial(appConfig.FileStore)) e.GET("/modeltype", handlers.ModelType(appConfig.FileStore)) e.GET("/modelversion", handlers.ModelVersion(appConfig.FileStore)) - e.GET("/geospatialdata", handlers.GeospatialData(appConfig)) + e.GET("/geospatialdata", handlers.GeospatialData(appConfig.FileStore, appConfig.DestinationCRS)) e.Logger.Fatal(e.Start(appConfig.Address())) } diff --git a/tools/geom.go b/tools/geom.go index 023eec2..7a2926f 100644 --- a/tools/geom.go +++ b/tools/geom.go @@ -2,17 +2,9 @@ package tools import ( "bufio" - "errors" - "fmt" - "math" "path/filepath" - "regexp" - "strconv" "strings" "sync" - - "github.com/USACE/filestore" - "github.com/dewberry/gdal" ) // GeomFileContents keywords and data container for ras flow file search @@ -25,293 +17,6 @@ type GeomFileContents struct { Structures []hydraulicStructures `json:"Hydraulic Structures"` } -type hydraulicStructures struct { - River string `json:"River Name"` - Reach string `json:"Reach Name"` - NumXS int `json:"Num CrossSections"` - NumCulverts int `json:"Num Culverts"` - BridgeData bridgeData `json:"Bridges"` - NumInlines int `json:"Num Inlines"` -} - -type bridgeData struct { - NumBridges int `json:"Num Bridges"` - Bridges []bridges `json:"Bridges"` -} - -type bridges struct { - Name string - Station float64 - Description string - DeckWidth float64 `json:"Deck Width"` - UpHighChord chordPairs `json:"Upstream High Chord"` - UpLowChord chordPairs `json:"Upstream Low Chord"` - DownHighChord chordPairs `json:"Downstream High Chord"` - DownLowChord chordPairs `json:"Downstream Max Chord"` - NumPiers int `json:"Num Piers"` -} - -type chordPairs struct { - Max float64 - Min float64 -} - -func maxValue(values []float64) (float64, error) { - if len(values) == 0 { - return 0.0, errors.New("Cannot detect a maximum value in an empty slice") - } - - max := values[0] - for _, v := range values { - if v > max { - max = v - } - } - - return max, nil -} - -func minValue(values []float64) (float64, error) { - if len(values) == 0 { - return 0.0, errors.New("Cannot detect a minimum value in an empty slice") - } - - min := values[0] - for _, v := range values { - if v < min { - min = v - } - } - - return min, nil -} - -func rightofEquals(line string) string { - return strings.TrimSpace(strings.Split(line, "=")[1]) -} - -func getDescription(sc *bufio.Scanner, idx int, endLine string) (string, int, error) { - description := "" - nLines := 0 - for sc.Scan() { - line := sc.Text() - idx++ - endDescription, err := regexp.MatchString(endLine, line) - if err != nil { - return description, idx, err - } - if endDescription { - return description, idx, nil - } - if line != "" { - if nLines > 0 { - description += "\n" - } - description += line - nLines++ - } - } - return description, idx, nil -} -func numberofLines(nValues int, colWidth int, valueWidth int) int { - nLines := math.Ceil(float64(nValues) / (float64(colWidth) / float64(valueWidth))) - return int(nLines) -} - -func datafromTextBlock(hsSc *bufio.Scanner, colWidth int, valueWidth int, nLines int, nSkipLines int) ([]float64, error) { - values := []float64{} - nSkipped := 0 - nProcessed := 0 -out: - for hsSc.Scan() { - if nSkipped < nSkipLines { - nSkipped++ - continue - } - nProcessed++ - line := hsSc.Text() - for s := 0; s < colWidth; { - if len(line) > s { - sVal := strings.TrimSpace(line[s : s+valueWidth]) - if sVal != "" { - val, err := strconv.ParseFloat(sVal, 64) - if err != nil { - return values, err - } - values = append(values, val) - } - s += valueWidth - } else { - if nLines == nProcessed { - break out - } - break - } - - } - if nLines == nProcessed { - break out - } - } - return values, nil -} - -func getHighLowChord(hsSc *bufio.Scanner, nElevsText string, colWidth int, valueWidth int) ([2]chordPairs, error) { - highLowPairs := [2]chordPairs{} - - nElevs, err := strconv.Atoi(strings.TrimSpace(nElevsText)) - if err != nil { - return highLowPairs, err - } - - nLines := numberofLines(nElevs, colWidth, valueWidth) - - elevHighChord, err := datafromTextBlock(hsSc, colWidth, valueWidth, nLines, nLines) - if err != nil { - return highLowPairs, err - } - - maxHighCord, err := maxValue(elevHighChord) - if err != nil { - return highLowPairs, err - } - - minHighCord, err := minValue(elevHighChord) - if err != nil { - return highLowPairs, err - } - highLowPairs[0] = chordPairs{Max: maxHighCord, Min: minHighCord} - - elevLowChord, err := datafromTextBlock(hsSc, 80, 8, nLines, 0) - if err != nil { - return highLowPairs, err - } - - maxLowCord, err := maxValue(elevLowChord) - if err != nil { - return highLowPairs, err - } - - minLowCord, err := minValue(elevLowChord) - if err != nil { - return highLowPairs, err - } - highLowPairs[1] = chordPairs{Max: maxLowCord, Min: minLowCord} - return highLowPairs, nil -} - -func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { - bridge := bridges{} - - station, err := strconv.ParseFloat(strings.TrimSpace(lineData[1]), 64) - if err != nil { - return bridge, err - } - bridge.Station = station - - for hsSc.Scan() { - line := hsSc.Text() - switch { - case strings.HasPrefix(line, "BEGIN DESCRIPTION"): - description, _, err := getDescription(hsSc, 0, "END DESCRIPTION:") - if err != nil { - return bridge, err - } - bridge.Description += description - - case strings.HasPrefix(line, "Node Name="): - bridge.Name = rightofEquals(line) - - case strings.HasPrefix(line, "Deck Dist"): - hsSc.Scan() - nextLineData := strings.Split(hsSc.Text(), ",") - deckWidth, err := strconv.ParseFloat(strings.TrimSpace(nextLineData[0]), 64) - if err != nil { - return bridge, err - } - bridge.DeckWidth = deckWidth - upHighLowPair, err := getHighLowChord(hsSc, nextLineData[4], 80, 8) - if err != nil { - return bridge, err - } - bridge.UpHighChord = upHighLowPair[0] - bridge.UpLowChord = upHighLowPair[1] - - downHighLowPair, err := getHighLowChord(hsSc, nextLineData[5], 80, 8) - if err != nil { - return bridge, err - } - bridge.DownHighChord = downHighLowPair[0] - bridge.DownLowChord = downHighLowPair[1] - - case strings.HasPrefix(line, "Pier Skew"): - bridge.NumPiers++ - - case strings.HasPrefix(line, "BR Coef"): - return bridge, err - } - } - return bridge, nil -} - -func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStructures, error) { - structures := hydraulicStructures{} - bData := bridgeData{} - - newf, err := rm.FileStore.GetObject(fn) - if err != nil { - return structures, nil - } - defer newf.Close() - - hsSc := bufio.NewScanner(newf) - - i := 0 - for hsSc.Scan() { - if i == idx { - riverReach := strings.Split(rightofEquals(hsSc.Text()), ",") - structures.River = strings.TrimSpace(riverReach[0]) - structures.Reach = strings.TrimSpace(riverReach[1]) - } else if i > idx { - line := hsSc.Text() - if strings.HasPrefix(line, "River Reach=") { - structures.BridgeData = bData - return structures, nil - } - if strings.HasPrefix(line, "Type RM Length L Ch R =") { - data := strings.Split(rightofEquals(line), ",") - structureType, err := strconv.Atoi(strings.TrimSpace(data[0])) - if err != nil { - return structures, err - } - switch structureType { - case 1: - structures.NumXS++ - - case 2: - structures.NumCulverts++ - - case 3: - bridge, err := getBridgeData(hsSc, data) - if err != nil { - return structures, err - } - bData.Bridges = append(bData.Bridges, bridge) - bData.NumBridges++ - - case 5: - structures.NumInlines++ - - } - } - } - i++ - } - structures.BridgeData = bData - - return structures, nil -} - // getGeomData Reads a geometry file. returns none to allow concurrency func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error) { @@ -368,414 +73,3 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error rm.Metadata.GeomFiles = append(rm.Metadata.GeomFiles, meta) return } - -// getProjection Reads a projection file. returns none to allow concurrency -func getProjection(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error) { - - defer wg.Done() - - f, err := rm.FileStore.GetObject(fn) - if err != nil { - errChan <- err - return - } - defer f.Close() - sc := bufio.NewScanner(f) - sc.Scan() - line := sc.Text() - if strings.HasPrefix(line, "PROJCS[") { - rm.Metadata.Projection = line - return - } - - errChan <- errors.New("Unexpected projection file structure") - return -} - -// GeoData ... -type GeoData struct { - Features map[string]Features - Georeference int -} - -// Features ... -type Features struct { - Rivers []vectorLayer - XS []vectorLayer - Banks []vectorLayer - StorageAreas []vectorLayer - TwoDAreas []vectorLayer - HydraulicStructures []vectorLayer -} - -type vectorLayer struct { - FeatureName string `json:"feature_name"` - Fields map[string]interface{} `json:"fields"` - Geometry []uint8 `json:"geometry"` -} - -type xyzPoint struct { - x float64 - y float64 - z float64 -} - -func dataPairsfromTextBlock(sc *bufio.Scanner, nPairs int, colWidth int, valueWidth int) ([][2]float64, error) { - var stride int = valueWidth * 2 - pairs := [][2]float64{} -out: - for sc.Scan() { - line := sc.Text() - for s := 0; s < colWidth; { - if len(line) > s { - val1, err := strconv.ParseFloat(strings.TrimSpace(line[s:s+valueWidth]), 64) - if err != nil { - return pairs, err - } - val2, err := strconv.ParseFloat(strings.TrimSpace(line[s+valueWidth:s+stride]), 64) - if err != nil { - return pairs, err - } - pairs = append(pairs, [2]float64{val1, val2}) - if len(pairs) == nPairs { - break out - } - } else { - break - } - s += stride - } - } - return pairs, nil -} - -func getDataPairsfromTextBlock(nDataPairsLine string, sc *bufio.Scanner, colWidth int, valueWidth int) ([][2]float64, error) { - pairs := [][2]float64{} - for sc.Scan() { - line := sc.Text() - if strings.HasPrefix(line, nDataPairsLine) { - nPairs, err := strconv.Atoi(rightofEquals(line)) - if err != nil { - return pairs, err - } - pairs, err = dataPairsfromTextBlock(sc, nPairs, colWidth, valueWidth) - if err != nil { - return pairs, err - } - break - } - } - return pairs, nil -} - -// distance returns the distance along a straight line in euclidean space -func distance(p0, p1 [2]float64) float64 { - result := math.Sqrt(math.Pow((p1[0]-p0[0]), 2) + math.Pow((p1[1]-p0[1]), 2)) - return result -} - -// pointAtDistance returns a new point along a straight line in euclidean space -// at a specified distance -func pointAtDistance(p0, p1 [2]float64, delta float64) [2]float64 { - distanceRatio := delta / distance(p0, p1) - newX := (1-distanceRatio)*p0[0] + distanceRatio*p1[0] - newY := (1-distanceRatio)*p0[1] + distanceRatio*p1[1] - return [2]float64{newX, newY} -} - -// interpZ creates a new point a given distance along a line composed -// of many segments. -func interpXY(xyPairs [][2]float64, d float64) [2]float64 { - // newPoint is an x, y pair - var newPoint [2]float64 - lineSegments := len(xyPairs) - 1 - lineLength := 0.0 - -findLineSegment: - for i := 0; i < lineSegments; i++ { - p0, p1 := xyPairs[i], xyPairs[i+1] - lineLength += distance(p0, p1) - - switch { - case lineLength > d: - delta := distance(p0, p1) - (lineLength - d) - newPoint = pointAtDistance(p0, p1, delta) - break findLineSegment - - default: - continue - } - } - if d >= lineLength { - if d-lineLength <= 0.1 { - return xyPairs[len(xyPairs)-1] - } - fmt.Printf("The interpolated point has a station of %v while the xy line is %v long", d, lineLength) - } - return newPoint -} - -// attributeZ using station from cross-section line and gis coordinates -func attributeZ(xyPairs [][2]float64, mzPairs [][2]float64) []xyzPoint { - points := []xyzPoint{} - startingStation := mzPairs[0][0] - for _, mzPair := range mzPairs { - newPoint := interpXY(xyPairs, mzPair[0]-startingStation) - if newPoint[0] != 0 && newPoint[1] != 0 { - points = append(points, xyzPoint{newPoint[0], newPoint[1], mzPair[1]}) - } else { - fmt.Printf("Interpolated point has an xy value of (%v, %v). ", newPoint[0], newPoint[1]) - } - } - return points -} - -func getTransform(sourceCRS string, destinationCRS int) (gdal.CoordinateTransform, error) { - transform := gdal.CoordinateTransform{} - sourceSpRef := gdal.CreateSpatialReference(sourceCRS) - sourceSpRef.MorphFromESRI() - if err := sourceSpRef.Validate(); err != nil { - return transform, errors.New("Unable to extract geospatial data. Invalid source Projection") - } - - destinationSpRef := gdal.CreateSpatialReference("") - if err := destinationSpRef.FromEPSG(destinationCRS); err != nil { - return transform, err - } - transform = gdal.CreateCoordinateTransform(sourceSpRef, destinationSpRef) - return transform, nil -} - -func flipXYLineString(xyLineString gdal.Geometry) gdal.Geometry { - yxLineString := gdal.Create(gdal.GT_LineString) - nPoints := xyLineString.PointCount() - for i := 0; i < nPoints; i++ { - x, y, _ := xyLineString.Point(i) - yxLineString.AddPoint2D(y, x) - } - xyLineString.Destroy() - return yxLineString -} - -func flipXYLineString25D(xyzLineString gdal.Geometry) gdal.Geometry { - yxzLineString := gdal.Create(gdal.GT_LineString25D) - nPoints := xyzLineString.PointCount() - for i := 0; i < nPoints; i++ { - x, y, z := xyzLineString.Point(i) - yxzLineString.AddPoint(y, x, z) - } - xyzLineString.Destroy() - return yxzLineString -} - -func flipXYLinearRing(xyLinearRing gdal.Geometry) gdal.Geometry { - yxLinearRing := gdal.Create(gdal.GT_LinearRing) - nPoints := xyLinearRing.PointCount() - for i := 0; i < nPoints; i++ { - x, y, _ := xyLinearRing.Point(i) - yxLinearRing.AddPoint2D(y, x) - } - xyLinearRing.Destroy() - return yxLinearRing -} - -func flipXYPoint(xyPoint gdal.Geometry) gdal.Geometry { - yxPoint := gdal.Create(gdal.GT_Point) - nPoints := xyPoint.PointCount() - for i := 0; i < nPoints; i++ { - x, y, _ := xyPoint.Point(i) - yxPoint.AddPoint2D(y, x) - } - xyPoint.Destroy() - return yxPoint -} - -func getRiverCenterline(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vectorLayer, error) { - riverReach := strings.Split(rightofEquals(sc.Text()), ",") - layer := vectorLayer{FeatureName: fmt.Sprintf("%s, %s", strings.TrimSpace(riverReach[0]), strings.TrimSpace(riverReach[1]))} - - xyPairs, err := getDataPairsfromTextBlock("Reach XY=", sc, 64, 16) - if err != nil { - return layer, err - } - - xyLineString := gdal.Create(gdal.GT_LineString) - for _, pair := range xyPairs { - xyLineString.AddPoint2D(pair[0], pair[1]) - } - - xyLineString.Transform(transform) - // This is a temporary fix since the x and y values need to be flipped: - yxLineString := flipXYLineString(xyLineString) - - multiLineString := yxLineString.ForceToMultiLineString() - wkb, err := multiLineString.ToWKB() - if err != nil { - return layer, err - } - layer.Geometry = wkb - return layer, err -} - -func getXSBanks(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName string) (vectorLayer, []vectorLayer, error) { - bankLayers := []vectorLayer{} - - xsLayer, xyPairs, startingStation, err := getXS(sc, transform, riverReachName) - if err != nil { - return xsLayer, bankLayers, err - } - for sc.Scan() { - line := sc.Text() - if strings.HasPrefix(line, "Bank Sta=") { - bankLayers, err = getBanks(line, transform, xsLayer, xyPairs, startingStation) - if err != nil { - return xsLayer, bankLayers, err - } - break - } - } - return xsLayer, bankLayers, err -} - -func getXS(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName string) (vectorLayer, [][2]float64, float64, error) { - compData := strings.Split(rightofEquals(sc.Text()), ",") - layer := vectorLayer{FeatureName: strings.TrimSpace(compData[1]), Fields: map[string]interface{}{}} - layer.Fields["RiverReachName"] = riverReachName - - xyPairs, err := getDataPairsfromTextBlock("XS GIS Cut Line", sc, 64, 16) - if err != nil { - return layer, xyPairs, 0.0, err - } - - mzPairs, err := getDataPairsfromTextBlock("#Sta/Elev", sc, 80, 8) - if err != nil { - return layer, xyPairs, mzPairs[0][0], err - } - - xyzPoints := attributeZ(xyPairs, mzPairs) - - xyzLineString := gdal.Create(gdal.GT_LineString25D) - for _, point := range xyzPoints { - xyzLineString.AddPoint(point.x, point.y, point.z) - } - - xyzLineString.Transform(transform) - // This is a temporary fix since the x and y values need to be flipped - yxzLineString := flipXYLineString25D(xyzLineString) - - multiLineString := yxzLineString.ForceToMultiLineString() - wkb, err := multiLineString.ToWKB() - if err != nil { - return layer, xyPairs, mzPairs[0][0], err - } - layer.Geometry = wkb - return layer, xyPairs, mzPairs[0][0], err -} - -func getBanks(line string, transform gdal.CoordinateTransform, xsLayer vectorLayer, xyPairs [][2]float64, startingStation float64) ([]vectorLayer, error) { - layers := []vectorLayer{} - - bankStations := strings.Split(rightofEquals(line), ",") - for _, s := range bankStations { - layer := vectorLayer{FeatureName: strings.TrimSpace(s), Fields: map[string]interface{}{}} - layer.Fields["RiverReachName"] = xsLayer.Fields["RiverReachName"] - layer.Fields["xsName"] = xsLayer.FeatureName - bankStation, err := strconv.ParseFloat(s, 64) - if err != nil { - return layers, err - } - bankXY := interpXY(xyPairs, bankStation-startingStation) - xyPoint := gdal.Create(gdal.GT_Point) - xyPoint.AddPoint2D(bankXY[0], bankXY[1]) - xyPoint.Transform(transform) - // This is a temporary fix since the x and y values need to be flipped - yxPoint := flipXYPoint(xyPoint) - multiPoint := yxPoint.ForceToMultiPoint() - wkb, err := multiPoint.ToWKB() - if err != nil { - return layers, err - } - layer.Geometry = wkb - layers = append(layers, layer) - } - return layers, nil -} - -func getStorageArea(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vectorLayer, error) { - layer := vectorLayer{FeatureName: strings.TrimSpace(strings.Split(rightofEquals(sc.Text()), ",")[0])} - - xyPairs, err := getDataPairsfromTextBlock("Storage Area Surface Line=", sc, 32, 16) - if err != nil { - return layer, err - } - - xyLinearRing := gdal.Create(gdal.GT_LinearRing) - for _, pair := range xyPairs { - xyLinearRing.AddPoint2D(pair[0], pair[1]) - } - - xyLinearRing.Transform(transform) - // This is a temporary fix since the x and y values need to be flipped: - yxLinearRing := flipXYLinearRing(xyLinearRing) - - yxPolygon := gdal.Create(gdal.GT_Polygon) - yxPolygon.AddGeometry(yxLinearRing) - yxMultiPolygon := yxPolygon.ForceToMultiPolygon() - wkb, err := yxMultiPolygon.ToWKB() - if err != nil { - return layer, err - } - layer.Geometry = wkb - return layer, err -} - -// GetGeospatialData ... -func GetGeospatialData(gd *GeoData, fs filestore.FileStore, geomFilePath string, sourceCRS string, destinationCRS int) error { - geomFileName := filepath.Base(geomFilePath) - f := Features{} - riverReachName := "" - - file, err := fs.GetObject(geomFilePath) - if err != nil { - return err - } - defer file.Close() - - sc := bufio.NewScanner(file) - - transform, err := getTransform(sourceCRS, destinationCRS) - if err != nil { - return err - } - - for sc.Scan() { - line := sc.Text() - switch { - case strings.HasPrefix(line, "River Reach="): - riverLayer, err := getRiverCenterline(sc, transform) - if err != nil { - return err - } - f.Rivers = append(f.Rivers, riverLayer) - riverReachName = riverLayer.FeatureName - - case strings.HasPrefix(line, "Storage Area="): - storageAreaLayer, err := getStorageArea(sc, transform) - if err != nil { - return err - } - f.StorageAreas = append(f.StorageAreas, storageAreaLayer) - - case strings.HasPrefix(line, "Type RM Length L Ch R = 1"): - xsLayer, bankLayers, err := getXSBanks(sc, transform, riverReachName) - if err != nil { - return err - } - f.XS = append(f.XS, xsLayer) - f.Banks = append(f.Banks, bankLayers...) - } - } - - gd.Features[geomFileName] = f - return nil -} diff --git a/tools/geospatial.go b/tools/geospatial.go new file mode 100644 index 0000000..34da024 --- /dev/null +++ b/tools/geospatial.go @@ -0,0 +1,401 @@ +package tools + +import ( + "bufio" + "errors" + "fmt" + "math" + "path/filepath" + "strconv" + "strings" + + "github.com/USACE/filestore" + "github.com/dewberry/gdal" +) + +// GeoData ... +type GeoData struct { + Features map[string]features + Georeference int +} + +type features struct { + Rivers []vectorLayer + XS []vectorLayer + Banks []vectorLayer + StorageAreas []vectorLayer + TwoDAreas []vectorLayer + HydraulicStructures []vectorLayer +} + +type vectorLayer struct { + FeatureName string `json:"feature_name"` + Fields map[string]interface{} `json:"fields"` + Geometry []uint8 `json:"geometry"` +} + +type xyzPoint struct { + x float64 + y float64 + z float64 +} + +func dataPairsfromTextBlock(sc *bufio.Scanner, nPairs int, colWidth int, valueWidth int) ([][2]float64, error) { + var stride int = valueWidth * 2 + pairs := [][2]float64{} +out: + for sc.Scan() { + line := sc.Text() + for s := 0; s < colWidth; { + if len(line) > s { + val1, err := strconv.ParseFloat(strings.TrimSpace(line[s:s+valueWidth]), 64) + if err != nil { + return pairs, err + } + val2, err := strconv.ParseFloat(strings.TrimSpace(line[s+valueWidth:s+stride]), 64) + if err != nil { + return pairs, err + } + pairs = append(pairs, [2]float64{val1, val2}) + if len(pairs) == nPairs { + break out + } + } else { + break + } + s += stride + } + } + return pairs, nil +} + +func getDataPairsfromTextBlock(nDataPairsLine string, sc *bufio.Scanner, colWidth int, valueWidth int) ([][2]float64, error) { + pairs := [][2]float64{} + for sc.Scan() { + line := sc.Text() + if strings.HasPrefix(line, nDataPairsLine) { + nPairs, err := strconv.Atoi(rightofEquals(line)) + if err != nil { + return pairs, err + } + pairs, err = dataPairsfromTextBlock(sc, nPairs, colWidth, valueWidth) + if err != nil { + return pairs, err + } + break + } + } + return pairs, nil +} + +// distance returns the distance along a straight line in euclidean space +func distance(p0, p1 [2]float64) float64 { + result := math.Sqrt(math.Pow((p1[0]-p0[0]), 2) + math.Pow((p1[1]-p0[1]), 2)) + return result +} + +// pointAtDistance returns a new point along a straight line in euclidean space +// at a specified distance +func pointAtDistance(p0, p1 [2]float64, delta float64) [2]float64 { + distanceRatio := delta / distance(p0, p1) + newX := (1-distanceRatio)*p0[0] + distanceRatio*p1[0] + newY := (1-distanceRatio)*p0[1] + distanceRatio*p1[1] + return [2]float64{newX, newY} +} + +// interpZ creates a new point a given distance along a line composed +// of many segments. +func interpXY(xyPairs [][2]float64, d float64) [2]float64 { + // newPoint is an x, y pair + var newPoint [2]float64 + lineSegments := len(xyPairs) - 1 + lineLength := 0.0 + +findLineSegment: + for i := 0; i < lineSegments; i++ { + p0, p1 := xyPairs[i], xyPairs[i+1] + lineLength += distance(p0, p1) + + switch { + case lineLength > d: + delta := distance(p0, p1) - (lineLength - d) + newPoint = pointAtDistance(p0, p1, delta) + break findLineSegment + + default: + continue + } + } + if d >= lineLength { + if d-lineLength <= 0.1 { + return xyPairs[len(xyPairs)-1] + } + fmt.Printf("The interpolated point has a station of %v while the xy line is %v long", d, lineLength) + } + return newPoint +} + +// attributeZ using station from cross-section line and gis coordinates +func attributeZ(xyPairs [][2]float64, mzPairs [][2]float64) []xyzPoint { + points := []xyzPoint{} + startingStation := mzPairs[0][0] + for _, mzPair := range mzPairs { + newPoint := interpXY(xyPairs, mzPair[0]-startingStation) + if newPoint[0] != 0 && newPoint[1] != 0 { + points = append(points, xyzPoint{newPoint[0], newPoint[1], mzPair[1]}) + } else { + fmt.Printf("Interpolated point has an xy value of (%v, %v). ", newPoint[0], newPoint[1]) + } + } + return points +} + +func getTransform(sourceCRS string, destinationCRS int) (gdal.CoordinateTransform, error) { + transform := gdal.CoordinateTransform{} + sourceSpRef := gdal.CreateSpatialReference(sourceCRS) + sourceSpRef.MorphFromESRI() + if err := sourceSpRef.Validate(); err != nil { + return transform, errors.New("Unable to extract geospatial data. Invalid source Projection") + } + + destinationSpRef := gdal.CreateSpatialReference("") + if err := destinationSpRef.FromEPSG(destinationCRS); err != nil { + return transform, err + } + transform = gdal.CreateCoordinateTransform(sourceSpRef, destinationSpRef) + return transform, nil +} + +func flipXYLineString(xyLineString gdal.Geometry) gdal.Geometry { + yxLineString := gdal.Create(gdal.GT_LineString) + nPoints := xyLineString.PointCount() + for i := 0; i < nPoints; i++ { + x, y, _ := xyLineString.Point(i) + yxLineString.AddPoint2D(y, x) + } + xyLineString.Destroy() + return yxLineString +} + +func flipXYLineString25D(xyzLineString gdal.Geometry) gdal.Geometry { + yxzLineString := gdal.Create(gdal.GT_LineString25D) + nPoints := xyzLineString.PointCount() + for i := 0; i < nPoints; i++ { + x, y, z := xyzLineString.Point(i) + yxzLineString.AddPoint(y, x, z) + } + xyzLineString.Destroy() + return yxzLineString +} + +func flipXYLinearRing(xyLinearRing gdal.Geometry) gdal.Geometry { + yxLinearRing := gdal.Create(gdal.GT_LinearRing) + nPoints := xyLinearRing.PointCount() + for i := 0; i < nPoints; i++ { + x, y, _ := xyLinearRing.Point(i) + yxLinearRing.AddPoint2D(y, x) + } + xyLinearRing.Destroy() + return yxLinearRing +} + +func flipXYPoint(xyPoint gdal.Geometry) gdal.Geometry { + yxPoint := gdal.Create(gdal.GT_Point) + nPoints := xyPoint.PointCount() + for i := 0; i < nPoints; i++ { + x, y, _ := xyPoint.Point(i) + yxPoint.AddPoint2D(y, x) + } + xyPoint.Destroy() + return yxPoint +} + +func getRiverCenterline(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vectorLayer, error) { + riverReach := strings.Split(rightofEquals(sc.Text()), ",") + layer := vectorLayer{FeatureName: fmt.Sprintf("%s, %s", strings.TrimSpace(riverReach[0]), strings.TrimSpace(riverReach[1]))} + + xyPairs, err := getDataPairsfromTextBlock("Reach XY=", sc, 64, 16) + if err != nil { + return layer, err + } + + xyLineString := gdal.Create(gdal.GT_LineString) + for _, pair := range xyPairs { + xyLineString.AddPoint2D(pair[0], pair[1]) + } + + xyLineString.Transform(transform) + // This is a temporary fix since the x and y values need to be flipped: + yxLineString := flipXYLineString(xyLineString) + + multiLineString := yxLineString.ForceToMultiLineString() + wkb, err := multiLineString.ToWKB() + if err != nil { + return layer, err + } + layer.Geometry = wkb + return layer, err +} + +func getXSBanks(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName string) (vectorLayer, []vectorLayer, error) { + bankLayers := []vectorLayer{} + + xsLayer, xyPairs, startingStation, err := getXS(sc, transform, riverReachName) + if err != nil { + return xsLayer, bankLayers, err + } + for sc.Scan() { + line := sc.Text() + if strings.HasPrefix(line, "Bank Sta=") { + bankLayers, err = getBanks(line, transform, xsLayer, xyPairs, startingStation) + if err != nil { + return xsLayer, bankLayers, err + } + break + } + } + return xsLayer, bankLayers, err +} + +func getXS(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName string) (vectorLayer, [][2]float64, float64, error) { + compData := strings.Split(rightofEquals(sc.Text()), ",") + layer := vectorLayer{FeatureName: strings.TrimSpace(compData[1]), Fields: map[string]interface{}{}} + layer.Fields["RiverReachName"] = riverReachName + + xyPairs, err := getDataPairsfromTextBlock("XS GIS Cut Line", sc, 64, 16) + if err != nil { + return layer, xyPairs, 0.0, err + } + + mzPairs, err := getDataPairsfromTextBlock("#Sta/Elev", sc, 80, 8) + if err != nil { + return layer, xyPairs, mzPairs[0][0], err + } + + xyzPoints := attributeZ(xyPairs, mzPairs) + + xyzLineString := gdal.Create(gdal.GT_LineString25D) + for _, point := range xyzPoints { + xyzLineString.AddPoint(point.x, point.y, point.z) + } + + xyzLineString.Transform(transform) + // This is a temporary fix since the x and y values need to be flipped + yxzLineString := flipXYLineString25D(xyzLineString) + + multiLineString := yxzLineString.ForceToMultiLineString() + wkb, err := multiLineString.ToWKB() + if err != nil { + return layer, xyPairs, mzPairs[0][0], err + } + layer.Geometry = wkb + return layer, xyPairs, mzPairs[0][0], err +} + +func getBanks(line string, transform gdal.CoordinateTransform, xsLayer vectorLayer, xyPairs [][2]float64, startingStation float64) ([]vectorLayer, error) { + layers := []vectorLayer{} + + bankStations := strings.Split(rightofEquals(line), ",") + for _, s := range bankStations { + layer := vectorLayer{FeatureName: strings.TrimSpace(s), Fields: map[string]interface{}{}} + layer.Fields["RiverReachName"] = xsLayer.Fields["RiverReachName"] + layer.Fields["xsName"] = xsLayer.FeatureName + bankStation, err := strconv.ParseFloat(s, 64) + if err != nil { + return layers, err + } + bankXY := interpXY(xyPairs, bankStation-startingStation) + xyPoint := gdal.Create(gdal.GT_Point) + xyPoint.AddPoint2D(bankXY[0], bankXY[1]) + xyPoint.Transform(transform) + // This is a temporary fix since the x and y values need to be flipped + yxPoint := flipXYPoint(xyPoint) + multiPoint := yxPoint.ForceToMultiPoint() + wkb, err := multiPoint.ToWKB() + if err != nil { + return layers, err + } + layer.Geometry = wkb + layers = append(layers, layer) + } + return layers, nil +} + +func getStorageArea(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vectorLayer, error) { + layer := vectorLayer{FeatureName: strings.TrimSpace(strings.Split(rightofEquals(sc.Text()), ",")[0])} + + xyPairs, err := getDataPairsfromTextBlock("Storage Area Surface Line=", sc, 32, 16) + if err != nil { + return layer, err + } + + xyLinearRing := gdal.Create(gdal.GT_LinearRing) + for _, pair := range xyPairs { + xyLinearRing.AddPoint2D(pair[0], pair[1]) + } + + xyLinearRing.Transform(transform) + // This is a temporary fix since the x and y values need to be flipped: + yxLinearRing := flipXYLinearRing(xyLinearRing) + + yxPolygon := gdal.Create(gdal.GT_Polygon) + yxPolygon.AddGeometry(yxLinearRing) + yxMultiPolygon := yxPolygon.ForceToMultiPolygon() + wkb, err := yxMultiPolygon.ToWKB() + if err != nil { + return layer, err + } + layer.Geometry = wkb + return layer, err +} + +// GetGeospatialData ... +func GetGeospatialData(gd *GeoData, fs filestore.FileStore, geomFilePath string, sourceCRS string, destinationCRS int) error { + geomFileName := filepath.Base(geomFilePath) + f := features{} + riverReachName := "" + + file, err := fs.GetObject(geomFilePath) + if err != nil { + return err + } + defer file.Close() + + sc := bufio.NewScanner(file) + + transform, err := getTransform(sourceCRS, destinationCRS) + if err != nil { + return err + } + + for sc.Scan() { + line := sc.Text() + switch { + case strings.HasPrefix(line, "River Reach="): + riverLayer, err := getRiverCenterline(sc, transform) + if err != nil { + return err + } + f.Rivers = append(f.Rivers, riverLayer) + riverReachName = riverLayer.FeatureName + + case strings.HasPrefix(line, "Storage Area="): + storageAreaLayer, err := getStorageArea(sc, transform) + if err != nil { + return err + } + f.StorageAreas = append(f.StorageAreas, storageAreaLayer) + + case strings.HasPrefix(line, "Type RM Length L Ch R = 1"): + xsLayer, bankLayers, err := getXSBanks(sc, transform, riverReachName) + if err != nil { + return err + } + f.XS = append(f.XS, xsLayer) + f.Banks = append(f.Banks, bankLayers...) + } + } + + gd.Features[geomFileName] = f + return nil +} diff --git a/tools/model.go b/tools/model.go index 63a9954..0e6b342 100644 --- a/tools/model.go +++ b/tools/model.go @@ -1,10 +1,12 @@ package tools import ( + "bufio" "errors" "fmt" "path/filepath" "regexp" + "strings" "sync" "github.com/USACE/filestore" @@ -54,6 +56,7 @@ type Model struct { Files ModelFiles } +// ModelFiles ... type ModelFiles struct { InputFiles InputFiles OutputFiles OutputFiles @@ -103,6 +106,7 @@ type SupplementalFiles struct { ObservationalData interface{} // placeholder } +// RasModel ... type RasModel struct { FileStore filestore.FileStore ModelDirectory string @@ -113,10 +117,12 @@ type RasModel struct { FileList []string } +// IsAModel ... func (rm *RasModel) IsAModel() bool { return rm.isModel } +// IsGeospatial ... func (rm *RasModel) IsGeospatial() bool { if rm.Metadata.GeomFiles[0].FileExt != "" { return true @@ -124,14 +130,17 @@ func (rm *RasModel) IsGeospatial() bool { return false } +// ModelType ... func (rm *RasModel) ModelType() string { return rm.Type } +// ModelVersion ... func (rm *RasModel) ModelVersion() string { return rm.Version } +// Index ... func (rm *RasModel) Index() (Model, error) { if !rm.IsAModel() { return Model{}, errors.New("model is not valid") @@ -193,7 +202,7 @@ func (rm *RasModel) Index() (Model, error) { // GeospatialData ... func (rm *RasModel) GeospatialData(destinationCRS int) (GeoData, error) { - gd := GeoData{Features: make(map[string]Features), Georeference: destinationCRS} + gd := GeoData{Features: make(map[string]features), Georeference: destinationCRS} sourceCRS := rm.Metadata.Projection for _, g := range rm.Metadata.GeomFiles { @@ -221,6 +230,30 @@ func getModelFiles(rm *RasModel) error { return nil } +// getProjection Reads a projection file. returns none to allow concurrency +func getProjection(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error) { + + defer wg.Done() + + f, err := rm.FileStore.GetObject(fn) + if err != nil { + errChan <- err + return + } + defer f.Close() + sc := bufio.NewScanner(f) + sc.Scan() + line := sc.Text() + if strings.HasPrefix(line, "PROJCS[") { + rm.Metadata.Projection = line + return + } + + errChan <- errors.New("Unexpected projection file structure") + return +} + +// NewRasModel ... func NewRasModel(key string, fs filestore.FileStore) (*RasModel, error) { rm := RasModel{ModelDirectory: filepath.Dir(key), FileStore: fs, Type: "RAS"} diff --git a/tools/structures.go b/tools/structures.go new file mode 100644 index 0000000..c557a3f --- /dev/null +++ b/tools/structures.go @@ -0,0 +1,238 @@ +package tools + +import ( + "bufio" + "math" + "strconv" + "strings" +) + +type hydraulicStructures struct { + River string `json:"River Name"` + Reach string `json:"Reach Name"` + NumXS int `json:"Num CrossSections"` + NumCulverts int `json:"Num Culverts"` + BridgeData bridgeData `json:"Bridges"` + NumInlines int `json:"Num Inlines"` +} + +type bridgeData struct { + NumBridges int `json:"Num Bridges"` + Bridges []bridges `json:"Bridges"` +} + +type bridges struct { + Name string + Station float64 + Description string + DeckWidth float64 `json:"Deck Width"` + UpHighChord chordPairs `json:"Upstream High Chord"` + UpLowChord chordPairs `json:"Upstream Low Chord"` + DownHighChord chordPairs `json:"Downstream High Chord"` + DownLowChord chordPairs `json:"Downstream Max Chord"` + NumPiers int `json:"Num Piers"` +} + +type chordPairs struct { + Max float64 + Min float64 +} + +func numberofLines(nValues int, colWidth int, valueWidth int) int { + nLines := math.Ceil(float64(nValues) / (float64(colWidth) / float64(valueWidth))) + return int(nLines) +} + +func datafromTextBlock(hsSc *bufio.Scanner, colWidth int, valueWidth int, nLines int, nSkipLines int) ([]float64, error) { + values := []float64{} + nSkipped := 0 + nProcessed := 0 +out: + for hsSc.Scan() { + if nSkipped < nSkipLines { + nSkipped++ + continue + } + nProcessed++ + line := hsSc.Text() + for s := 0; s < colWidth; { + if len(line) > s { + sVal := strings.TrimSpace(line[s : s+valueWidth]) + if sVal != "" { + val, err := strconv.ParseFloat(sVal, 64) + if err != nil { + return values, err + } + values = append(values, val) + } + s += valueWidth + } else { + if nLines == nProcessed { + break out + } + break + } + + } + if nLines == nProcessed { + break out + } + } + return values, nil +} + +func getHighLowChord(hsSc *bufio.Scanner, nElevsText string, colWidth int, valueWidth int) ([2]chordPairs, error) { + highLowPairs := [2]chordPairs{} + + nElevs, err := strconv.Atoi(strings.TrimSpace(nElevsText)) + if err != nil { + return highLowPairs, err + } + + nLines := numberofLines(nElevs, colWidth, valueWidth) + + elevHighChord, err := datafromTextBlock(hsSc, colWidth, valueWidth, nLines, nLines) + if err != nil { + return highLowPairs, err + } + + maxHighCord, err := maxValue(elevHighChord) + if err != nil { + return highLowPairs, err + } + + minHighCord, err := minValue(elevHighChord) + if err != nil { + return highLowPairs, err + } + highLowPairs[0] = chordPairs{Max: maxHighCord, Min: minHighCord} + + elevLowChord, err := datafromTextBlock(hsSc, 80, 8, nLines, 0) + if err != nil { + return highLowPairs, err + } + + maxLowCord, err := maxValue(elevLowChord) + if err != nil { + return highLowPairs, err + } + + minLowCord, err := minValue(elevLowChord) + if err != nil { + return highLowPairs, err + } + highLowPairs[1] = chordPairs{Max: maxLowCord, Min: minLowCord} + return highLowPairs, nil +} + +func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { + bridge := bridges{} + + station, err := strconv.ParseFloat(strings.TrimSpace(lineData[1]), 64) + if err != nil { + return bridge, err + } + bridge.Station = station + + for hsSc.Scan() { + line := hsSc.Text() + switch { + case strings.HasPrefix(line, "BEGIN DESCRIPTION"): + description, _, err := getDescription(hsSc, 0, "END DESCRIPTION:") + if err != nil { + return bridge, err + } + bridge.Description += description + + case strings.HasPrefix(line, "Node Name="): + bridge.Name = rightofEquals(line) + + case strings.HasPrefix(line, "Deck Dist"): + hsSc.Scan() + nextLineData := strings.Split(hsSc.Text(), ",") + deckWidth, err := strconv.ParseFloat(strings.TrimSpace(nextLineData[0]), 64) + if err != nil { + return bridge, err + } + bridge.DeckWidth = deckWidth + upHighLowPair, err := getHighLowChord(hsSc, nextLineData[4], 80, 8) + if err != nil { + return bridge, err + } + bridge.UpHighChord = upHighLowPair[0] + bridge.UpLowChord = upHighLowPair[1] + + downHighLowPair, err := getHighLowChord(hsSc, nextLineData[5], 80, 8) + if err != nil { + return bridge, err + } + bridge.DownHighChord = downHighLowPair[0] + bridge.DownLowChord = downHighLowPair[1] + + case strings.HasPrefix(line, "Pier Skew"): + bridge.NumPiers++ + + case strings.HasPrefix(line, "BR Coef"): + return bridge, err + } + } + return bridge, nil +} + +func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStructures, error) { + structures := hydraulicStructures{} + bData := bridgeData{} + + newf, err := rm.FileStore.GetObject(fn) + if err != nil { + return structures, nil + } + defer newf.Close() + + hsSc := bufio.NewScanner(newf) + + i := 0 + for hsSc.Scan() { + if i == idx { + riverReach := strings.Split(rightofEquals(hsSc.Text()), ",") + structures.River = strings.TrimSpace(riverReach[0]) + structures.Reach = strings.TrimSpace(riverReach[1]) + } else if i > idx { + line := hsSc.Text() + if strings.HasPrefix(line, "River Reach=") { + structures.BridgeData = bData + return structures, nil + } + if strings.HasPrefix(line, "Type RM Length L Ch R =") { + data := strings.Split(rightofEquals(line), ",") + structureType, err := strconv.Atoi(strings.TrimSpace(data[0])) + if err != nil { + return structures, err + } + switch structureType { + case 1: + structures.NumXS++ + + case 2: + structures.NumCulverts++ + + case 3: + bridge, err := getBridgeData(hsSc, data) + if err != nil { + return structures, err + } + bData.Bridges = append(bData.Bridges, bridge) + bData.NumBridges++ + + case 5: + structures.NumInlines++ + + } + } + } + i++ + } + structures.BridgeData = bData + + return structures, nil +} diff --git a/tools/utils.go b/tools/utils.go new file mode 100644 index 0000000..4273d82 --- /dev/null +++ b/tools/utils.go @@ -0,0 +1,66 @@ +package tools + +import ( + "bufio" + "errors" + "regexp" + "strings" +) + +func maxValue(values []float64) (float64, error) { + if len(values) == 0 { + return 0.0, errors.New("Cannot detect a maximum value in an empty slice") + } + + max := values[0] + for _, v := range values { + if v > max { + max = v + } + } + + return max, nil +} + +func minValue(values []float64) (float64, error) { + if len(values) == 0 { + return 0.0, errors.New("Cannot detect a minimum value in an empty slice") + } + + min := values[0] + for _, v := range values { + if v < min { + min = v + } + } + + return min, nil +} + +func rightofEquals(line string) string { + return strings.TrimSpace(strings.Split(line, "=")[1]) +} + +func getDescription(sc *bufio.Scanner, idx int, endLine string) (string, int, error) { + description := "" + nLines := 0 + for sc.Scan() { + line := sc.Text() + idx++ + endDescription, err := regexp.MatchString(endLine, line) + if err != nil { + return description, idx, err + } + if endDescription { + return description, idx, nil + } + if line != "" { + if nLines > 0 { + description += "\n" + } + description += line + nLines++ + } + } + return description, idx, nil +} From 2e5ff22c93719ac6009052d03d53d54eead077b2 Mon Sep 17 00:00:00 2001 From: Shane Date: Mon, 4 Jan 2021 15:19:21 -0500 Subject: [PATCH 08/14] Export Features and VectorLayer types --- tools/geospatial.go | 44 +++++++++++++++++++++++--------------------- tools/model.go | 2 +- 2 files changed, 24 insertions(+), 22 deletions(-) diff --git a/tools/geospatial.go b/tools/geospatial.go index 34da024..17daca9 100644 --- a/tools/geospatial.go +++ b/tools/geospatial.go @@ -15,20 +15,22 @@ import ( // GeoData ... type GeoData struct { - Features map[string]features + Features map[string]Features Georeference int } -type features struct { - Rivers []vectorLayer - XS []vectorLayer - Banks []vectorLayer - StorageAreas []vectorLayer - TwoDAreas []vectorLayer - HydraulicStructures []vectorLayer +// Features ... +type Features struct { + Rivers []VectorLayer + XS []VectorLayer + Banks []VectorLayer + StorageAreas []VectorLayer + TwoDAreas []VectorLayer + HydraulicStructures []VectorLayer } -type vectorLayer struct { +// VectorLayer ... +type VectorLayer struct { FeatureName string `json:"feature_name"` Fields map[string]interface{} `json:"fields"` Geometry []uint8 `json:"geometry"` @@ -210,9 +212,9 @@ func flipXYPoint(xyPoint gdal.Geometry) gdal.Geometry { return yxPoint } -func getRiverCenterline(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vectorLayer, error) { +func getRiverCenterline(sc *bufio.Scanner, transform gdal.CoordinateTransform) (VectorLayer, error) { riverReach := strings.Split(rightofEquals(sc.Text()), ",") - layer := vectorLayer{FeatureName: fmt.Sprintf("%s, %s", strings.TrimSpace(riverReach[0]), strings.TrimSpace(riverReach[1]))} + layer := VectorLayer{FeatureName: fmt.Sprintf("%s, %s", strings.TrimSpace(riverReach[0]), strings.TrimSpace(riverReach[1]))} xyPairs, err := getDataPairsfromTextBlock("Reach XY=", sc, 64, 16) if err != nil { @@ -237,8 +239,8 @@ func getRiverCenterline(sc *bufio.Scanner, transform gdal.CoordinateTransform) ( return layer, err } -func getXSBanks(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName string) (vectorLayer, []vectorLayer, error) { - bankLayers := []vectorLayer{} +func getXSBanks(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName string) (VectorLayer, []VectorLayer, error) { + bankLayers := []VectorLayer{} xsLayer, xyPairs, startingStation, err := getXS(sc, transform, riverReachName) if err != nil { @@ -257,9 +259,9 @@ func getXSBanks(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReac return xsLayer, bankLayers, err } -func getXS(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName string) (vectorLayer, [][2]float64, float64, error) { +func getXS(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName string) (VectorLayer, [][2]float64, float64, error) { compData := strings.Split(rightofEquals(sc.Text()), ",") - layer := vectorLayer{FeatureName: strings.TrimSpace(compData[1]), Fields: map[string]interface{}{}} + layer := VectorLayer{FeatureName: strings.TrimSpace(compData[1]), Fields: map[string]interface{}{}} layer.Fields["RiverReachName"] = riverReachName xyPairs, err := getDataPairsfromTextBlock("XS GIS Cut Line", sc, 64, 16) @@ -292,12 +294,12 @@ func getXS(sc *bufio.Scanner, transform gdal.CoordinateTransform, riverReachName return layer, xyPairs, mzPairs[0][0], err } -func getBanks(line string, transform gdal.CoordinateTransform, xsLayer vectorLayer, xyPairs [][2]float64, startingStation float64) ([]vectorLayer, error) { - layers := []vectorLayer{} +func getBanks(line string, transform gdal.CoordinateTransform, xsLayer VectorLayer, xyPairs [][2]float64, startingStation float64) ([]VectorLayer, error) { + layers := []VectorLayer{} bankStations := strings.Split(rightofEquals(line), ",") for _, s := range bankStations { - layer := vectorLayer{FeatureName: strings.TrimSpace(s), Fields: map[string]interface{}{}} + layer := VectorLayer{FeatureName: strings.TrimSpace(s), Fields: map[string]interface{}{}} layer.Fields["RiverReachName"] = xsLayer.Fields["RiverReachName"] layer.Fields["xsName"] = xsLayer.FeatureName bankStation, err := strconv.ParseFloat(s, 64) @@ -321,8 +323,8 @@ func getBanks(line string, transform gdal.CoordinateTransform, xsLayer vectorLay return layers, nil } -func getStorageArea(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vectorLayer, error) { - layer := vectorLayer{FeatureName: strings.TrimSpace(strings.Split(rightofEquals(sc.Text()), ",")[0])} +func getStorageArea(sc *bufio.Scanner, transform gdal.CoordinateTransform) (VectorLayer, error) { + layer := VectorLayer{FeatureName: strings.TrimSpace(strings.Split(rightofEquals(sc.Text()), ",")[0])} xyPairs, err := getDataPairsfromTextBlock("Storage Area Surface Line=", sc, 32, 16) if err != nil { @@ -352,7 +354,7 @@ func getStorageArea(sc *bufio.Scanner, transform gdal.CoordinateTransform) (vect // GetGeospatialData ... func GetGeospatialData(gd *GeoData, fs filestore.FileStore, geomFilePath string, sourceCRS string, destinationCRS int) error { geomFileName := filepath.Base(geomFilePath) - f := features{} + f := Features{} riverReachName := "" file, err := fs.GetObject(geomFilePath) diff --git a/tools/model.go b/tools/model.go index 0e6b342..f1edc1c 100644 --- a/tools/model.go +++ b/tools/model.go @@ -202,7 +202,7 @@ func (rm *RasModel) Index() (Model, error) { // GeospatialData ... func (rm *RasModel) GeospatialData(destinationCRS int) (GeoData, error) { - gd := GeoData{Features: make(map[string]features), Georeference: destinationCRS} + gd := GeoData{Features: make(map[string]Features), Georeference: destinationCRS} sourceCRS := rm.Metadata.Projection for _, g := range rm.Metadata.GeomFiles { From 324fcd2f722e3932cc92c898e50e8df96ef5a9d2 Mon Sep 17 00:00:00 2001 From: Shane Date: Tue, 5 Jan 2021 14:29:11 -0500 Subject: [PATCH 09/14] Extract culvert structures --- tools/structures.go | 236 +++++++++++++++++++++++++++++++++++++++----- tools/utils.go | 7 +- 2 files changed, 212 insertions(+), 31 deletions(-) diff --git a/tools/structures.go b/tools/structures.go index c557a3f..9f7a218 100644 --- a/tools/structures.go +++ b/tools/structures.go @@ -7,13 +7,52 @@ import ( "strings" ) +var conduitShapes map[int]string = map[int]string{ + 1: "Circular", + 2: "Box", + 3: "Pipe Arch", + 4: "Ellipse", + 5: "Arch", + 6: "Semi-Circle", + 7: "Low Arch", + 8: "High Arch", + 9: "Conspan Arch"} + type hydraulicStructures struct { - River string `json:"River Name"` - Reach string `json:"Reach Name"` - NumXS int `json:"Num CrossSections"` + River string `json:"River Name"` + Reach string `json:"Reach Name"` + NumXS int `json:"Num CrossSections"` + CulvertData culvertData `json:"Culverts"` + BridgeData bridgeData `json:"Bridges"` + NumInlines int `json:"Num Inlines"` +} + +type culvertData struct { NumCulverts int `json:"Num Culverts"` - BridgeData bridgeData `json:"Bridges"` - NumInlines int `json:"Num Inlines"` + Culverts []culverts `json:"Culverts"` +} + +type culverts struct { + Name string + Station float64 + Description string + DeckWidth float64 `json:"Deck Width"` + UpHighChord chordPairs `json:"Upstream High Chord"` + UpLowChord chordPairs `json:"Upstream Low Chord"` + DownHighChord chordPairs `json:"Downstream High Chord"` + DownLowChord chordPairs `json:"Downstream Low Chord"` + NumConduits int `json:"Num Culvert Conduits"` + Conduits []conduits `json:"Culvert Conduits"` +} + +type conduits struct { + NumBarrels int `json:"Num Barrels"` + Group string `json:"Culvert Group"` + Shape string + Rise float64 + Span float64 + Length float64 + ManningsN float64 `json:"Mannings N"` } type bridgeData struct { @@ -29,7 +68,7 @@ type bridges struct { UpHighChord chordPairs `json:"Upstream High Chord"` UpLowChord chordPairs `json:"Upstream Low Chord"` DownHighChord chordPairs `json:"Downstream High Chord"` - DownLowChord chordPairs `json:"Downstream Max Chord"` + DownLowChord chordPairs `json:"Downstream Low Chord"` NumPiers int `json:"Num Piers"` } @@ -43,7 +82,7 @@ func numberofLines(nValues int, colWidth int, valueWidth int) int { return int(nLines) } -func datafromTextBlock(hsSc *bufio.Scanner, colWidth int, valueWidth int, nLines int, nSkipLines int) ([]float64, error) { +func datafromTextBlock(hsSc *bufio.Scanner, nLines int, nSkipLines int, colWidth int, valueWidth int) ([]float64, error) { values := []float64{} nSkipped := 0 nProcessed := 0 @@ -81,48 +120,186 @@ out: return values, nil } -func getHighLowChord(hsSc *bufio.Scanner, nElevsText string, colWidth int, valueWidth int) ([2]chordPairs, error) { +func getChord(hsSc *bufio.Scanner, nLines int, nSkipLines int, colWidth int, valueWidth int) (chordPairs, error) { + pair := chordPairs{} + + elevChord, err := datafromTextBlock(hsSc, nLines, nSkipLines, colWidth, valueWidth) + + if err != nil { + return pair, err + } + + if len(elevChord) == 0 { + return pair, nil + } + + maxChord, err := maxValue(elevChord) + if err != nil { + return pair, err + } + + minChord, err := minValue(elevChord) + if err != nil { + return pair, err + } + + pair = chordPairs{Max: maxChord, Min: minChord} + return pair, nil +} + +func getHighLowChord(hsSc *bufio.Scanner, nElevText string, colWidth int, valueWidth int) ([2]chordPairs, error) { highLowPairs := [2]chordPairs{} - nElevs, err := strconv.Atoi(strings.TrimSpace(nElevsText)) + nElev, err := strconv.Atoi(strings.TrimSpace(nElevText)) if err != nil { return highLowPairs, err } - nLines := numberofLines(nElevs, colWidth, valueWidth) + nLines := numberofLines(nElev, 80, 8) - elevHighChord, err := datafromTextBlock(hsSc, colWidth, valueWidth, nLines, nLines) + highPair, err := getChord(hsSc, nLines, nLines, 80, 8) if err != nil { return highLowPairs, err } + highLowPairs[0] = highPair - maxHighCord, err := maxValue(elevHighChord) + lowPair, err := getChord(hsSc, nLines, 0, 80, 8) if err != nil { return highLowPairs, err } + highLowPairs[1] = lowPair + + return highLowPairs, nil +} + +func stringtoFloat(s string) (float64, error) { + trimmed := strings.TrimSpace(s) + if trimmed != "" { + val, err := strconv.ParseFloat(trimmed, 64) + if err != nil { + return 0, err + } + return val, nil + } + return 0, nil +} - minHighCord, err := minValue(elevHighChord) +func getConduits(line string, single bool, shapeMap map[int]string) (conduits, error) { + lineData := strings.Split(rightofEquals(line), ",") + conduit := conduits{} + + if single { + conduit.NumBarrels = 1 + conduit.Group = strings.TrimSpace(lineData[13]) + + } else { + numbarrels, err := strconv.Atoi(strings.TrimSpace(lineData[11])) + if err != nil { + return conduit, err + } + conduit.NumBarrels = numbarrels + conduit.Group = strings.TrimSpace(lineData[12]) + } + + shapeID, err := strconv.Atoi(strings.TrimSpace(lineData[0])) if err != nil { - return highLowPairs, err + return conduit, err } - highLowPairs[0] = chordPairs{Max: maxHighCord, Min: minHighCord} + conduit.Shape = shapeMap[shapeID] - elevLowChord, err := datafromTextBlock(hsSc, 80, 8, nLines, 0) + rise, err := stringtoFloat(lineData[1]) if err != nil { - return highLowPairs, err + return conduit, err } + conduit.Rise = rise - maxLowCord, err := maxValue(elevLowChord) + span, err := stringtoFloat(lineData[2]) if err != nil { - return highLowPairs, err + return conduit, err } + conduit.Span = span - minLowCord, err := minValue(elevLowChord) + length, err := stringtoFloat(lineData[3]) if err != nil { - return highLowPairs, err + return conduit, err } - highLowPairs[1] = chordPairs{Max: maxLowCord, Min: minLowCord} - return highLowPairs, nil + conduit.Length = length + + mannings, err := stringtoFloat(lineData[4]) + if err != nil { + return conduit, err + } + conduit.ManningsN = mannings + + return conduit, nil +} + +func getCulvertData(hsSc *bufio.Scanner, lineData []string, shapeMap map[int]string) (culverts, error) { + culvert := culverts{} + + station, err := strconv.ParseFloat(strings.TrimSpace(lineData[1]), 64) + if err != nil { + return culvert, err + } + culvert.Station = station + + for hsSc.Scan() { + line := hsSc.Text() + switch { + case strings.HasPrefix(line, "BEGIN DESCRIPTION"): + description, _, err := getDescription(hsSc, 0, "END DESCRIPTION:") + if err != nil { + return culvert, err + } + culvert.Description += description + + case strings.HasPrefix(line, "Node Name="): + culvert.Name = rightofEquals(line) + + case strings.HasPrefix(line, "Deck Dist"): + hsSc.Scan() + nextLineData := strings.Split(hsSc.Text(), ",") + deckWidth, err := strconv.ParseFloat(strings.TrimSpace(nextLineData[0]), 64) + if err != nil { + return culvert, err + } + culvert.DeckWidth = deckWidth + + upHighLowPair, err := getHighLowChord(hsSc, nextLineData[4], 80, 8) + if err != nil { + return culvert, err + } + culvert.UpHighChord = upHighLowPair[0] + culvert.UpLowChord = upHighLowPair[1] + + downHighLowPair, err := getHighLowChord(hsSc, nextLineData[5], 80, 8) + if err != nil { + return culvert, err + } + culvert.DownHighChord = downHighLowPair[0] + culvert.DownLowChord = downHighLowPair[1] + + case strings.HasPrefix(line, "Culvert="): + conduit, err := getConduits(line, true, shapeMap) + if err != nil { + return culvert, err + } + culvert.Conduits = append(culvert.Conduits, conduit) + culvert.NumConduits++ + + case strings.HasPrefix(line, "Multiple Barrel Culv="): + conduit, err := getConduits(line, false, shapeMap) + if err != nil { + return culvert, err + } + culvert.Conduits = append(culvert.Conduits, conduit) + culvert.NumConduits++ + + case strings.HasPrefix(line, "BC Design"): + return culvert, err + } + } + return culvert, nil } func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { @@ -155,6 +332,7 @@ func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { return bridge, err } bridge.DeckWidth = deckWidth + upHighLowPair, err := getHighLowChord(hsSc, nextLineData[4], 80, 8) if err != nil { return bridge, err @@ -172,7 +350,7 @@ func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { case strings.HasPrefix(line, "Pier Skew"): bridge.NumPiers++ - case strings.HasPrefix(line, "BR Coef"): + case strings.HasPrefix(line, "BC Design"): return bridge, err } } @@ -182,6 +360,7 @@ func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStructures, error) { structures := hydraulicStructures{} bData := bridgeData{} + cData := culvertData{} newf, err := rm.FileStore.GetObject(fn) if err != nil { @@ -200,6 +379,7 @@ func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStruc } else if i > idx { line := hsSc.Text() if strings.HasPrefix(line, "River Reach=") { + structures.CulvertData = cData structures.BridgeData = bData return structures, nil } @@ -214,7 +394,12 @@ func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStruc structures.NumXS++ case 2: - structures.NumCulverts++ + culvert, err := getCulvertData(hsSc, data, conduitShapes) + if err != nil { + return structures, err + } + cData.Culverts = append(cData.Culverts, culvert) + cData.NumCulverts++ case 3: bridge, err := getBridgeData(hsSc, data) @@ -232,6 +417,7 @@ func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStruc } i++ } + structures.CulvertData = cData structures.BridgeData = bData return structures, nil diff --git a/tools/utils.go b/tools/utils.go index 4273d82..7ff086b 100644 --- a/tools/utils.go +++ b/tools/utils.go @@ -3,7 +3,6 @@ package tools import ( "bufio" "errors" - "regexp" "strings" ) @@ -47,11 +46,7 @@ func getDescription(sc *bufio.Scanner, idx int, endLine string) (string, int, er for sc.Scan() { line := sc.Text() idx++ - endDescription, err := regexp.MatchString(endLine, line) - if err != nil { - return description, idx, err - } - if endDescription { + if strings.HasPrefix(line, endLine) { return description, idx, nil } if line != "" { From e313132e0c050b4f6d857d36984292bdae812bed Mon Sep 17 00:00:00 2001 From: Shane Date: Wed, 6 Jan 2021 17:17:25 -0500 Subject: [PATCH 10/14] Extract inline weirs --- tools/geom.go | 7 +- tools/structures.go | 353 +++++++++++++++++++++++++++++++++----------- tools/utils.go | 2 +- 3 files changed, 268 insertions(+), 94 deletions(-) diff --git a/tools/geom.go b/tools/geom.go index 7a2926f..9865173 100644 --- a/tools/geom.go +++ b/tools/geom.go @@ -2,6 +2,7 @@ package tools import ( "bufio" + "fmt" "path/filepath" "strings" "sync" @@ -38,6 +39,7 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error header := true idx := 0 for sc.Scan() { + idx++ line := sc.Text() switch { case strings.HasPrefix(line, "Geom Title="): @@ -59,8 +61,8 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error case strings.HasPrefix(line, "River Reach="): structures, err := getHydraulicStructureData(rm, fn, idx) if err != nil { - errChan <- err - return + fmt.Println(err) + continue } meta.Structures = append(meta.Structures, structures) header = false @@ -68,7 +70,6 @@ func getGeomData(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan error case strings.HasPrefix(line, "Storage Area="): header = false } - idx++ } rm.Metadata.GeomFiles = append(rm.Metadata.GeomFiles, meta) return diff --git a/tools/structures.go b/tools/structures.go index 9f7a218..0a52c4a 100644 --- a/tools/structures.go +++ b/tools/structures.go @@ -2,6 +2,7 @@ package tools import ( "bufio" + "errors" "math" "strconv" "strings" @@ -22,9 +23,9 @@ type hydraulicStructures struct { River string `json:"River Name"` Reach string `json:"Reach Name"` NumXS int `json:"Num CrossSections"` - CulvertData culvertData `json:"Culverts"` - BridgeData bridgeData `json:"Bridges"` - NumInlines int `json:"Num Inlines"` + CulvertData culvertData `json:"Culvert Data"` + BridgeData bridgeData `json:"Bridge Data"` + WeirData weirData `json:"Inline Weir Data"` } type culvertData struct { @@ -36,18 +37,23 @@ type culverts struct { Name string Station float64 Description string - DeckWidth float64 `json:"Deck Width"` - UpHighChord chordPairs `json:"Upstream High Chord"` - UpLowChord chordPairs `json:"Upstream Low Chord"` - DownHighChord chordPairs `json:"Downstream High Chord"` - DownLowChord chordPairs `json:"Downstream Low Chord"` - NumConduits int `json:"Num Culvert Conduits"` - Conduits []conduits `json:"Culvert Conduits"` + DeckWidth float64 `json:"Deck Width"` + UpHighChord maxMinPairs `json:"Upstream High Chord"` + UpLowChord maxMinPairs `json:"Upstream Low Chord"` + DownHighChord maxMinPairs `json:"Downstream High Chord"` + DownLowChord maxMinPairs `json:"Downstream Low Chord"` + NumConduits int `json:"Num Culvert Conduits"` + Conduits []conduits `json:"Culvert Conduits"` +} + +type maxMinPairs struct { + Max float64 + Min float64 } type conduits struct { - NumBarrels int `json:"Num Barrels"` - Group string `json:"Culvert Group"` + Name string + NumBarrels int `json:"Num Barrels"` Shape string Rise float64 Span float64 @@ -64,30 +70,46 @@ type bridges struct { Name string Station float64 Description string - DeckWidth float64 `json:"Deck Width"` - UpHighChord chordPairs `json:"Upstream High Chord"` - UpLowChord chordPairs `json:"Upstream Low Chord"` - DownHighChord chordPairs `json:"Downstream High Chord"` - DownLowChord chordPairs `json:"Downstream Low Chord"` - NumPiers int `json:"Num Piers"` + DeckWidth float64 `json:"Deck Width"` + UpHighChord maxMinPairs `json:"Upstream High Chord"` + UpLowChord maxMinPairs `json:"Upstream Low Chord"` + DownHighChord maxMinPairs `json:"Downstream High Chord"` + DownLowChord maxMinPairs `json:"Downstream Low Chord"` + NumPiers int `json:"Num Piers"` } -type chordPairs struct { - Max float64 - Min float64 +type weirData struct { + NumWeirs int `json:"Num Inline Weirs"` + Weirs []weirs `json:"Inline Weirs"` } -func numberofLines(nValues int, colWidth int, valueWidth int) int { - nLines := math.Ceil(float64(nValues) / (float64(colWidth) / float64(valueWidth))) - return int(nLines) +type weirs struct { + Name string + Station float64 + Description string + WeirWidth float64 `json:"Weir Width"` + WeirElev maxMinPairs `json:"Weir Elevations"` + NumGates int `json:"Num Gates"` + Gates []gates + NumConduits int `json:"Num Culvert Conduits"` + Conduits []conduits `json:"Culvert Conduits"` +} + +type gates struct { + Name string + Width float64 + Height float64 + NumOpenings int `json:"Num Openings"` } -func datafromTextBlock(hsSc *bufio.Scanner, nLines int, nSkipLines int, colWidth int, valueWidth int) ([]float64, error) { +func datafromTextBlock(hsSc *bufio.Scanner, i int, nLines int, nSkipLines int, colWidth int, valueWidth int, interval int) ([]float64, int, error) { values := []float64{} nSkipped := 0 nProcessed := 0 + nvalues := 0 out: for hsSc.Scan() { + i++ if nSkipped < nSkipLines { nSkipped++ continue @@ -98,11 +120,14 @@ out: if len(line) > s { sVal := strings.TrimSpace(line[s : s+valueWidth]) if sVal != "" { - val, err := strconv.ParseFloat(sVal, 64) - if err != nil { - return values, err + nvalues++ + if nvalues%interval == 0 { + val, err := strconv.ParseFloat(sVal, 64) + if err != nil { + return values, i, err + } + values = append(values, val) } - values = append(values, val) } s += valueWidth } else { @@ -117,59 +142,64 @@ out: break out } } - return values, nil + return values, i, nil } -func getChord(hsSc *bufio.Scanner, nLines int, nSkipLines int, colWidth int, valueWidth int) (chordPairs, error) { - pair := chordPairs{} +func getMaxMinElev(hsSc *bufio.Scanner, i int, nLines int, nSkipLines int, colWidth int, valueWidth int, interval int) (maxMinPairs, int, error) { + pair := maxMinPairs{} - elevChord, err := datafromTextBlock(hsSc, nLines, nSkipLines, colWidth, valueWidth) + elevations, i, err := datafromTextBlock(hsSc, i, nLines, nSkipLines, colWidth, valueWidth, interval) if err != nil { - return pair, err + return pair, i, err } - if len(elevChord) == 0 { - return pair, nil + if len(elevations) == 0 { + return pair, i, nil } - maxChord, err := maxValue(elevChord) + maxElev, err := maxValue(elevations) if err != nil { - return pair, err + return pair, i, err } - minChord, err := minValue(elevChord) + minElev, err := minValue(elevations) if err != nil { - return pair, err + return pair, i, err } - pair = chordPairs{Max: maxChord, Min: minChord} - return pair, nil + pair = maxMinPairs{Max: maxElev, Min: minElev} + return pair, i, nil } -func getHighLowChord(hsSc *bufio.Scanner, nElevText string, colWidth int, valueWidth int) ([2]chordPairs, error) { - highLowPairs := [2]chordPairs{} +func numberofLines(nValues int, colWidth int, valueWidth int) int { + nLines := math.Ceil(float64(nValues) / (float64(colWidth) / float64(valueWidth))) + return int(nLines) +} + +func getHighLowChord(hsSc *bufio.Scanner, i int, nElevText string, colWidth int, valueWidth int) ([2]maxMinPairs, int, error) { + highLowPairs := [2]maxMinPairs{} nElev, err := strconv.Atoi(strings.TrimSpace(nElevText)) if err != nil { - return highLowPairs, err + return highLowPairs, i, err } nLines := numberofLines(nElev, 80, 8) - highPair, err := getChord(hsSc, nLines, nLines, 80, 8) + highPair, i, err := getMaxMinElev(hsSc, i, nLines, nLines, 80, 8, 1) if err != nil { - return highLowPairs, err + return highLowPairs, i, err } highLowPairs[0] = highPair - lowPair, err := getChord(hsSc, nLines, 0, 80, 8) + lowPair, i, err := getMaxMinElev(hsSc, i, nLines, 0, 80, 8, 1) if err != nil { - return highLowPairs, err + return highLowPairs, i, err } highLowPairs[1] = lowPair - return highLowPairs, nil + return highLowPairs, i, nil } func stringtoFloat(s string) (float64, error) { @@ -190,7 +220,7 @@ func getConduits(line string, single bool, shapeMap map[int]string) (conduits, e if single { conduit.NumBarrels = 1 - conduit.Group = strings.TrimSpace(lineData[13]) + conduit.Name = strings.TrimSpace(lineData[13]) } else { numbarrels, err := strconv.Atoi(strings.TrimSpace(lineData[11])) @@ -198,7 +228,7 @@ func getConduits(line string, single bool, shapeMap map[int]string) (conduits, e return conduit, err } conduit.NumBarrels = numbarrels - conduit.Group = strings.TrimSpace(lineData[12]) + conduit.Name = strings.TrimSpace(lineData[12]) } shapeID, err := strconv.Atoi(strings.TrimSpace(lineData[0])) @@ -234,22 +264,24 @@ func getConduits(line string, single bool, shapeMap map[int]string) (conduits, e return conduit, nil } -func getCulvertData(hsSc *bufio.Scanner, lineData []string, shapeMap map[int]string) (culverts, error) { +func getCulvertData(hsSc *bufio.Scanner, i int, lineData []string, shapeMap map[int]string) (culverts, int, error) { culvert := culverts{} station, err := strconv.ParseFloat(strings.TrimSpace(lineData[1]), 64) if err != nil { - return culvert, err + return culvert, i, err } culvert.Station = station for hsSc.Scan() { + i++ line := hsSc.Text() switch { case strings.HasPrefix(line, "BEGIN DESCRIPTION"): - description, _, err := getDescription(hsSc, 0, "END DESCRIPTION:") + var description string + description, i, err = getDescription(hsSc, i, "END DESCRIPTION:") if err != nil { - return culvert, err + return culvert, i, err } culvert.Description += description @@ -257,24 +289,27 @@ func getCulvertData(hsSc *bufio.Scanner, lineData []string, shapeMap map[int]str culvert.Name = rightofEquals(line) case strings.HasPrefix(line, "Deck Dist"): + i++ hsSc.Scan() nextLineData := strings.Split(hsSc.Text(), ",") deckWidth, err := strconv.ParseFloat(strings.TrimSpace(nextLineData[0]), 64) if err != nil { - return culvert, err + return culvert, i, err } culvert.DeckWidth = deckWidth - upHighLowPair, err := getHighLowChord(hsSc, nextLineData[4], 80, 8) + var upHighLowPair [2]maxMinPairs + upHighLowPair, i, err = getHighLowChord(hsSc, i, nextLineData[4], 80, 8) if err != nil { - return culvert, err + return culvert, i, err } culvert.UpHighChord = upHighLowPair[0] culvert.UpLowChord = upHighLowPair[1] - downHighLowPair, err := getHighLowChord(hsSc, nextLineData[5], 80, 8) + var downHighLowPair [2]maxMinPairs + downHighLowPair, i, err = getHighLowChord(hsSc, i, nextLineData[5], 80, 8) if err != nil { - return culvert, err + return culvert, i, err } culvert.DownHighChord = downHighLowPair[0] culvert.DownLowChord = downHighLowPair[1] @@ -282,7 +317,7 @@ func getCulvertData(hsSc *bufio.Scanner, lineData []string, shapeMap map[int]str case strings.HasPrefix(line, "Culvert="): conduit, err := getConduits(line, true, shapeMap) if err != nil { - return culvert, err + return culvert, i, err } culvert.Conduits = append(culvert.Conduits, conduit) culvert.NumConduits++ @@ -290,34 +325,42 @@ func getCulvertData(hsSc *bufio.Scanner, lineData []string, shapeMap map[int]str case strings.HasPrefix(line, "Multiple Barrel Culv="): conduit, err := getConduits(line, false, shapeMap) if err != nil { - return culvert, err + return culvert, i, err } culvert.Conduits = append(culvert.Conduits, conduit) culvert.NumConduits++ case strings.HasPrefix(line, "BC Design"): - return culvert, err + return culvert, i, nil + + case strings.HasPrefix(line, "Type RM Length L Ch R ="): + return culvert, i, errors.New("Failed to terminate parsing of culvert at 'BC Design'") + + case strings.HasPrefix(line, "River Reach="): + return culvert, i, errors.New("Failed to terminate parsing of culvert at 'BC Design'") } } - return culvert, nil + return culvert, i, nil } -func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { +func getBridgeData(hsSc *bufio.Scanner, i int, lineData []string) (bridges, int, error) { bridge := bridges{} station, err := strconv.ParseFloat(strings.TrimSpace(lineData[1]), 64) if err != nil { - return bridge, err + return bridge, i, err } bridge.Station = station for hsSc.Scan() { + i++ line := hsSc.Text() switch { case strings.HasPrefix(line, "BEGIN DESCRIPTION"): - description, _, err := getDescription(hsSc, 0, "END DESCRIPTION:") + var description string + description, i, err = getDescription(hsSc, i, "END DESCRIPTION:") if err != nil { - return bridge, err + return bridge, i, err } bridge.Description += description @@ -325,24 +368,27 @@ func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { bridge.Name = rightofEquals(line) case strings.HasPrefix(line, "Deck Dist"): + i++ hsSc.Scan() nextLineData := strings.Split(hsSc.Text(), ",") deckWidth, err := strconv.ParseFloat(strings.TrimSpace(nextLineData[0]), 64) if err != nil { - return bridge, err + return bridge, i, err } bridge.DeckWidth = deckWidth - upHighLowPair, err := getHighLowChord(hsSc, nextLineData[4], 80, 8) + var upHighLowPair [2]maxMinPairs + upHighLowPair, i, err = getHighLowChord(hsSc, i, nextLineData[4], 80, 8) if err != nil { - return bridge, err + return bridge, i, err } bridge.UpHighChord = upHighLowPair[0] bridge.UpLowChord = upHighLowPair[1] - downHighLowPair, err := getHighLowChord(hsSc, nextLineData[5], 80, 8) + var downHighLowPair [2]maxMinPairs + downHighLowPair, i, err = getHighLowChord(hsSc, i, nextLineData[5], 80, 8) if err != nil { - return bridge, err + return bridge, i, err } bridge.DownHighChord = downHighLowPair[0] bridge.DownLowChord = downHighLowPair[1] @@ -351,38 +397,153 @@ func getBridgeData(hsSc *bufio.Scanner, lineData []string) (bridges, error) { bridge.NumPiers++ case strings.HasPrefix(line, "BC Design"): - return bridge, err + return bridge, i, nil + + case strings.HasPrefix(line, "Type RM Length L Ch R ="): + return bridge, i, errors.New("Failed to terminate parsing of bridge at 'BC Design'") + + case strings.HasPrefix(line, "River Reach="): + return bridge, i, errors.New("Failed to terminate parsing of bridge at 'BC Design'") + } + } + return bridge, i, nil +} + +func getGates(nextLine string) (gates, error) { + gate := gates{} + + nextLineData := strings.Split(nextLine, ",") + + gate.Name = strings.TrimSpace(nextLineData[0]) + + width, err := stringtoFloat(nextLineData[1]) + if err != nil { + return gate, err + } + gate.Width = width + + height, err := stringtoFloat(nextLineData[2]) + if err != nil { + return gate, err + } + gate.Height = height + + numopenings, err := strconv.Atoi(strings.TrimSpace(nextLineData[13])) + if err != nil { + return gate, err + } + gate.NumOpenings = numopenings + + return gate, nil +} + +func getWeirData(rm *RasModel, fn string, i int, shapeMap map[int]string) (weirs, error) { + weir := weirs{} + + f, err := rm.FileStore.GetObject(fn) + if err != nil { + return weir, err + } + defer f.Close() + + wSc := bufio.NewScanner(f) + + wi := 0 + for wSc.Scan() { + wi++ + if wi == i { + lineData := strings.Split(rightofEquals(wSc.Text()), ",") + station, err := strconv.ParseFloat(strings.TrimSpace(lineData[1]), 64) + if err != nil { + return weir, err + } + weir.Station = station + } else if wi > i { + line := wSc.Text() + switch { + case strings.HasPrefix(line, "BEGIN DESCRIPTION"): + description, _, err := getDescription(wSc, 0, "END DESCRIPTION:") + if err != nil { + return weir, err + } + weir.Description += description + + case strings.HasPrefix(line, "Node Name="): + weir.Name = rightofEquals(line) + + case strings.HasPrefix(line, "#Inline Weir SE="): + nElev, err := strconv.Atoi(strings.TrimSpace(rightofEquals(line))) + if err != nil { + return weir, err + } + nLines := numberofLines(nElev*2, 80, 8) + + elev, _, err := getMaxMinElev(wSc, 0, nLines, 0, 80, 8, 2) + if err != nil { + return weir, err + } + weir.WeirElev = elev + + case strings.HasPrefix(line, "IW Dist,WD"): + wSc.Scan() + nextLineData := strings.Split(wSc.Text(), ",") + weirWidth, err := strconv.ParseFloat(strings.TrimSpace(nextLineData[1]), 64) + if err != nil { + return weir, err + } + weir.WeirWidth = weirWidth + + case strings.HasPrefix(line, "IW Gate Name"): + wSc.Scan() + gate, err := getGates(wSc.Text()) + if err != nil { + return weir, err + } + weir.Gates = append(weir.Gates, gate) + weir.NumGates++ + + case strings.HasPrefix(line, "IW Culv="): + conduit, err := getConduits(line, false, shapeMap) + if err != nil { + return weir, err + } + weir.Conduits = append(weir.Conduits, conduit) + weir.NumConduits++ + + case strings.HasPrefix(line, "Type RM Length L Ch R ="): + return weir, nil + + case strings.HasPrefix(line, "River Reach="): + return weir, nil + } } } - return bridge, nil + return weir, nil } func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStructures, error) { structures := hydraulicStructures{} bData := bridgeData{} cData := culvertData{} + wData := weirData{} - newf, err := rm.FileStore.GetObject(fn) + f, err := rm.FileStore.GetObject(fn) if err != nil { - return structures, nil + return structures, err } - defer newf.Close() + defer f.Close() - hsSc := bufio.NewScanner(newf) + hsSc := bufio.NewScanner(f) i := 0 for hsSc.Scan() { + i++ if i == idx { riverReach := strings.Split(rightofEquals(hsSc.Text()), ",") structures.River = strings.TrimSpace(riverReach[0]) structures.Reach = strings.TrimSpace(riverReach[1]) } else if i > idx { line := hsSc.Text() - if strings.HasPrefix(line, "River Reach=") { - structures.CulvertData = cData - structures.BridgeData = bData - return structures, nil - } if strings.HasPrefix(line, "Type RM Length L Ch R =") { data := strings.Split(rightofEquals(line), ",") structureType, err := strconv.Atoi(strings.TrimSpace(data[0])) @@ -394,7 +555,8 @@ func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStruc structures.NumXS++ case 2: - culvert, err := getCulvertData(hsSc, data, conduitShapes) + var culvert culverts + culvert, i, err = getCulvertData(hsSc, i, data, conduitShapes) if err != nil { return structures, err } @@ -402,7 +564,8 @@ func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStruc cData.NumCulverts++ case 3: - bridge, err := getBridgeData(hsSc, data) + var bridge bridges + bridge, i, err = getBridgeData(hsSc, i, data) if err != nil { return structures, err } @@ -410,15 +573,25 @@ func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStruc bData.NumBridges++ case 5: - structures.NumInlines++ - + weir, err := getWeirData(rm, fn, i, conduitShapes) + if err != nil { + return structures, err + } + wData.Weirs = append(wData.Weirs, weir) + wData.NumWeirs++ } } + if strings.HasPrefix(line, "River Reach=") { + structures.CulvertData = cData + structures.BridgeData = bData + structures.WeirData = wData + return structures, nil + } } - i++ } structures.CulvertData = cData structures.BridgeData = bData + structures.WeirData = wData return structures, nil } diff --git a/tools/utils.go b/tools/utils.go index 7ff086b..f404458 100644 --- a/tools/utils.go +++ b/tools/utils.go @@ -44,8 +44,8 @@ func getDescription(sc *bufio.Scanner, idx int, endLine string) (string, int, er description := "" nLines := 0 for sc.Scan() { - line := sc.Text() idx++ + line := sc.Text() if strings.HasPrefix(line, endLine) { return description, idx, nil } From 0562965f825a4a13aad8062be2ebcea48fde7263 Mon Sep 17 00:00:00 2001 From: Shane Date: Thu, 7 Jan 2021 12:11:24 -0500 Subject: [PATCH 11/14] Check that the model and the coordinate reference system have the same units --- tools/geospatial.go | 21 +++++++++++++++++++++ tools/model.go | 17 ++++++++++++++++- 2 files changed, 37 insertions(+), 1 deletion(-) diff --git a/tools/geospatial.go b/tools/geospatial.go index 17daca9..4d832a8 100644 --- a/tools/geospatial.go +++ b/tools/geospatial.go @@ -42,6 +42,27 @@ type xyzPoint struct { z float64 } +var unitConsistencyMap map[string]string = map[string]string{ + "English Units": "US survey foot", + "SI Units": "metre"} + +// checkUnitConsistency checks that the unit system used by the model and its coordinate reference system are the same +func checkUnitConsistency(modelUnits string, sourceCRS string, ucMap map[string]string) error { + sourceSpRef := gdal.CreateSpatialReference(sourceCRS) + sourceSpRef.MorphFromESRI() + if err := sourceSpRef.Validate(); err != nil { + return errors.New("Unable to check unit consistency, invalid coordinate reference system") + } + + if crsUnits, ok := sourceSpRef.AttrValue("UNIT", 0); ok { + if ucMap[modelUnits] == crsUnits { + return nil + } + return errors.New("The unit system of the model and coordinate reference system are inconsistent") + } + return errors.New("Unable to check unit consistency, could not identify the coordinate reference system's units") +} + func dataPairsfromTextBlock(sc *bufio.Scanner, nPairs int, colWidth int, valueWidth int) ([][2]float64, error) { var stride int = valueWidth * 2 pairs := [][2]float64{} diff --git a/tools/model.go b/tools/model.go index f1edc1c..c512e59 100644 --- a/tools/model.go +++ b/tools/model.go @@ -202,15 +202,30 @@ func (rm *RasModel) Index() (Model, error) { // GeospatialData ... func (rm *RasModel) GeospatialData(destinationCRS int) (GeoData, error) { - gd := GeoData{Features: make(map[string]Features), Georeference: destinationCRS} + gd := GeoData{} + + modelUnits := rm.Metadata.ProjFileContents.Units sourceCRS := rm.Metadata.Projection + + if sourceCRS == "" { + return gd, errors.New("Cannot extract geospatial data, no coordinate reference system") + } + + if err := checkUnitConsistency(modelUnits, sourceCRS, unitConsistencyMap); err != nil { + return gd, err + } + + gd.Features = make(map[string]Features) + gd.Georeference = destinationCRS + for _, g := range rm.Metadata.GeomFiles { if err := GetGeospatialData(&gd, rm.FileStore, g.Path, sourceCRS, destinationCRS); err != nil { return gd, err } } return gd, nil + } func getModelFiles(rm *RasModel) error { From bac76fdcc271889c6c3e8736573093949173a6d1 Mon Sep 17 00:00:00 2001 From: Shane Date: Thu, 7 Jan 2021 16:35:53 -0500 Subject: [PATCH 12/14] Identify .prj or .projection as projection files --- tools/geospatial.go | 8 -------- tools/model.go | 33 ++++++++++++++++++--------------- tools/project.go | 3 --- 3 files changed, 18 insertions(+), 26 deletions(-) diff --git a/tools/geospatial.go b/tools/geospatial.go index 4d832a8..09c71bc 100644 --- a/tools/geospatial.go +++ b/tools/geospatial.go @@ -49,10 +49,6 @@ var unitConsistencyMap map[string]string = map[string]string{ // checkUnitConsistency checks that the unit system used by the model and its coordinate reference system are the same func checkUnitConsistency(modelUnits string, sourceCRS string, ucMap map[string]string) error { sourceSpRef := gdal.CreateSpatialReference(sourceCRS) - sourceSpRef.MorphFromESRI() - if err := sourceSpRef.Validate(); err != nil { - return errors.New("Unable to check unit consistency, invalid coordinate reference system") - } if crsUnits, ok := sourceSpRef.AttrValue("UNIT", 0); ok { if ucMap[modelUnits] == crsUnits { @@ -176,10 +172,6 @@ func attributeZ(xyPairs [][2]float64, mzPairs [][2]float64) []xyzPoint { func getTransform(sourceCRS string, destinationCRS int) (gdal.CoordinateTransform, error) { transform := gdal.CoordinateTransform{} sourceSpRef := gdal.CreateSpatialReference(sourceCRS) - sourceSpRef.MorphFromESRI() - if err := sourceSpRef.Validate(); err != nil { - return transform, errors.New("Unable to extract geospatial data. Invalid source Projection") - } destinationSpRef := gdal.CreateSpatialReference("") if err := destinationSpRef.FromEPSG(destinationCRS); err != nil { diff --git a/tools/model.go b/tools/model.go index c512e59..cd5aebf 100644 --- a/tools/model.go +++ b/tools/model.go @@ -6,10 +6,10 @@ import ( "fmt" "path/filepath" "regexp" - "strings" "sync" "github.com/USACE/filestore" + "github.com/dewberry/gdal" ) type fileExtMatchers struct { @@ -37,7 +37,7 @@ var rasRE fileExtMatchers = fileExtMatchers{ // Maybe these ones are better? nee SteadyRun: regexp.MustCompile(".r[0-9][0-9]"), // `^\.r(0[1-9]|[1-9][0-9])$` UnsteadyRun: regexp.MustCompile(".x[0-9][0-9]"), // `^\.x(0[1-9]|[1-9][0-9])$` AllFlowRun: regexp.MustCompile(".[rx][0-9][0-9]"), // `^\.[rx](0[1-9]|[1-9][0-9])$` - Projection: regexp.MustCompile(".projection"), + Projection: regexp.MustCompile(".pr[oj]"), } // holder of multiple wait groups to help process files concurrency @@ -209,7 +209,7 @@ func (rm *RasModel) GeospatialData(destinationCRS int) (GeoData, error) { sourceCRS := rm.Metadata.Projection if sourceCRS == "" { - return gd, errors.New("Cannot extract geospatial data, no coordinate reference system") + return gd, errors.New("Cannot extract geospatial data, no valid coordinate reference system") } if err := checkUnitConsistency(modelUnits, sourceCRS, unitConsistencyMap); err != nil { @@ -237,9 +237,7 @@ func getModelFiles(rm *RasModel) error { } for _, file := range *files { - if filepath.Ext(file.Name) != ".prj" { - rm.FileList = append(rm.FileList, filepath.Join(file.Path, file.Name)) - } + rm.FileList = append(rm.FileList, filepath.Join(file.Path, file.Name)) } return nil @@ -256,15 +254,23 @@ func getProjection(rm *RasModel, fn string, wg *sync.WaitGroup, errChan chan err return } defer f.Close() + sc := bufio.NewScanner(f) sc.Scan() line := sc.Text() - if strings.HasPrefix(line, "PROJCS[") { - rm.Metadata.Projection = line + + sourceSpRef := gdal.CreateSpatialReference(line) + if err := sourceSpRef.Validate(); err != nil { + rm.Metadata.Projection = "" return } + if rm.Metadata.Projection != "" { + errChan <- errors.New("Multiple projection files identified, cannot determine coordinate reference system") + return + } + + rm.Metadata.Projection = line - errChan <- errors.New("Unexpected projection file structure") return } @@ -290,7 +296,6 @@ func NewRasModel(key string, fs filestore.FileStore) (*RasModel, error) { errChan := make(chan error) var rasWG rasWaitGroup - nProjection := 0 for _, fp := range rm.FileList { ext := filepath.Ext(fp) @@ -310,12 +315,10 @@ func NewRasModel(key string, fs filestore.FileStore) (*RasModel, error) { go getFlowData(&rm, fp, &rasWG.Flow, errChan) case rasRE.Projection.MatchString(ext): - if nProjection > 0 { - return &rm, errors.New("multiple projection files identified") + if filepath.Base(key) != filepath.Base(fp) { + rasWG.Projection.Add(1) + go getProjection(&rm, fp, &rasWG.Projection, errChan) } - rasWG.Projection.Add(1) - go getProjection(&rm, fp, &rasWG.Projection, errChan) - nProjection++ } } diff --git a/tools/project.go b/tools/project.go index d60402c..c890d36 100644 --- a/tools/project.go +++ b/tools/project.go @@ -53,7 +53,6 @@ func rmNewLineChar(s string) string { // verifyPrjPath identifies the .prj file within the passed model directory ... func verifyPrjPath(key string, rm *RasModel) error { - rasFiles := make([]string, 0) if filepath.Ext(key) != ".prj" { return fmt.Errorf("%s is not a .prj file", key) @@ -67,9 +66,7 @@ func verifyPrjPath(key string, rm *RasModel) error { return fmt.Errorf("%s is not a RAS Project file", key) } - rasFiles = append(rasFiles, key) rm.Metadata.ProjFilePath = key - rm.FileList = rasFiles return nil } From dd2c5b00f7f796fb15bf487fdcbe7397b4cb1d4e Mon Sep 17 00:00:00 2001 From: Shane Date: Mon, 11 Jan 2021 11:03:49 -0500 Subject: [PATCH 13/14] Remove empty version from RAS MCAT model struct --- tools/model.go | 15 ++++++++++++--- 1 file changed, 12 insertions(+), 3 deletions(-) diff --git a/tools/model.go b/tools/model.go index cd5aebf..04587d3 100644 --- a/tools/model.go +++ b/tools/model.go @@ -334,13 +334,22 @@ func NewRasModel(key string, fs filestore.FileStore) (*RasModel, error) { } for _, p := range rm.Metadata.PlanFiles { - rm.Version += fmt.Sprintf("%s: %s, ", p.FileExt, p.ProgramVersion) + version := p.ProgramVersion + if version != "" { + rm.Version += fmt.Sprintf("%s: %s, ", p.FileExt, version) + } } for _, g := range rm.Metadata.GeomFiles { - rm.Version += fmt.Sprintf("%s: %s, ", g.FileExt, g.ProgramVersion) + version := g.ProgramVersion + if version != "" { + rm.Version += fmt.Sprintf("%s: %s, ", g.FileExt, version) + } } for _, f := range rm.Metadata.FlowFiles { - rm.Version += fmt.Sprintf("%s: %s, ", f.FileExt, f.ProgramVersion) + version := f.ProgramVersion + if version != "" { + rm.Version += fmt.Sprintf("%s: %s, ", f.FileExt, version) + } } if len(rm.Version) >= 2 { From 99e9c67affb6ff381770a2a93adb4e12b221bd84 Mon Sep 17 00:00:00 2001 From: Shane Date: Mon, 11 Jan 2021 13:32:53 -0500 Subject: [PATCH 14/14] Pass appConfig through geospatialdata handler and avoid passing global variables --- handlers/geospatialdata.go | 8 ++++---- main.go | 2 +- tools/geospatial.go | 4 ++-- tools/model.go | 2 +- tools/structures.go | 18 +++++++++--------- 5 files changed, 17 insertions(+), 17 deletions(-) diff --git a/handlers/geospatialdata.go b/handlers/geospatialdata.go index 505e46b..3aeada3 100644 --- a/handlers/geospatialdata.go +++ b/handlers/geospatialdata.go @@ -3,7 +3,7 @@ package handlers import ( "net/http" - "github.com/USACE/filestore" + "github.com/USACE/mcat-ras/config" ras "github.com/USACE/mcat-ras/tools" "github.com/labstack/echo/v4" @@ -19,17 +19,17 @@ import ( // @Success 200 {object} interface{} // @Failure 500 {object} SimpleResponse // @Router /geospatialdata [get] -func GeospatialData(fs *filestore.FileStore, destinationCRS int) echo.HandlerFunc { +func GeospatialData(ac *config.APIConfig) echo.HandlerFunc { return func(c echo.Context) error { definitionFile := c.QueryParam("definition_file") - rm, err := ras.NewRasModel(definitionFile, *fs) + rm, err := ras.NewRasModel(definitionFile, *ac.FileStore) if err != nil { return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) } - data, err := rm.GeospatialData(destinationCRS) + data, err := rm.GeospatialData(ac.DestinationCRS) if err != nil { return c.JSON(http.StatusInternalServerError, SimpleResponse{http.StatusInternalServerError, err.Error()}) } diff --git a/main.go b/main.go index b6b0361..d413225 100644 --- a/main.go +++ b/main.go @@ -39,7 +39,7 @@ func main() { e.GET("/isgeospatial", handlers.IsGeospatial(appConfig.FileStore)) e.GET("/modeltype", handlers.ModelType(appConfig.FileStore)) e.GET("/modelversion", handlers.ModelVersion(appConfig.FileStore)) - e.GET("/geospatialdata", handlers.GeospatialData(appConfig.FileStore, appConfig.DestinationCRS)) + e.GET("/geospatialdata", handlers.GeospatialData(appConfig)) e.Logger.Fatal(e.Start(appConfig.Address())) } diff --git a/tools/geospatial.go b/tools/geospatial.go index 09c71bc..143410e 100644 --- a/tools/geospatial.go +++ b/tools/geospatial.go @@ -47,11 +47,11 @@ var unitConsistencyMap map[string]string = map[string]string{ "SI Units": "metre"} // checkUnitConsistency checks that the unit system used by the model and its coordinate reference system are the same -func checkUnitConsistency(modelUnits string, sourceCRS string, ucMap map[string]string) error { +func checkUnitConsistency(modelUnits string, sourceCRS string) error { sourceSpRef := gdal.CreateSpatialReference(sourceCRS) if crsUnits, ok := sourceSpRef.AttrValue("UNIT", 0); ok { - if ucMap[modelUnits] == crsUnits { + if unitConsistencyMap[modelUnits] == crsUnits { return nil } return errors.New("The unit system of the model and coordinate reference system are inconsistent") diff --git a/tools/model.go b/tools/model.go index 04587d3..ac4d366 100644 --- a/tools/model.go +++ b/tools/model.go @@ -212,7 +212,7 @@ func (rm *RasModel) GeospatialData(destinationCRS int) (GeoData, error) { return gd, errors.New("Cannot extract geospatial data, no valid coordinate reference system") } - if err := checkUnitConsistency(modelUnits, sourceCRS, unitConsistencyMap); err != nil { + if err := checkUnitConsistency(modelUnits, sourceCRS); err != nil { return gd, err } diff --git a/tools/structures.go b/tools/structures.go index 0a52c4a..12fe2d8 100644 --- a/tools/structures.go +++ b/tools/structures.go @@ -214,7 +214,7 @@ func stringtoFloat(s string) (float64, error) { return 0, nil } -func getConduits(line string, single bool, shapeMap map[int]string) (conduits, error) { +func getConduits(line string, single bool) (conduits, error) { lineData := strings.Split(rightofEquals(line), ",") conduit := conduits{} @@ -235,7 +235,7 @@ func getConduits(line string, single bool, shapeMap map[int]string) (conduits, e if err != nil { return conduit, err } - conduit.Shape = shapeMap[shapeID] + conduit.Shape = conduitShapes[shapeID] rise, err := stringtoFloat(lineData[1]) if err != nil { @@ -264,7 +264,7 @@ func getConduits(line string, single bool, shapeMap map[int]string) (conduits, e return conduit, nil } -func getCulvertData(hsSc *bufio.Scanner, i int, lineData []string, shapeMap map[int]string) (culverts, int, error) { +func getCulvertData(hsSc *bufio.Scanner, i int, lineData []string) (culverts, int, error) { culvert := culverts{} station, err := strconv.ParseFloat(strings.TrimSpace(lineData[1]), 64) @@ -315,7 +315,7 @@ func getCulvertData(hsSc *bufio.Scanner, i int, lineData []string, shapeMap map[ culvert.DownLowChord = downHighLowPair[1] case strings.HasPrefix(line, "Culvert="): - conduit, err := getConduits(line, true, shapeMap) + conduit, err := getConduits(line, true) if err != nil { return culvert, i, err } @@ -323,7 +323,7 @@ func getCulvertData(hsSc *bufio.Scanner, i int, lineData []string, shapeMap map[ culvert.NumConduits++ case strings.HasPrefix(line, "Multiple Barrel Culv="): - conduit, err := getConduits(line, false, shapeMap) + conduit, err := getConduits(line, false) if err != nil { return culvert, i, err } @@ -437,7 +437,7 @@ func getGates(nextLine string) (gates, error) { return gate, nil } -func getWeirData(rm *RasModel, fn string, i int, shapeMap map[int]string) (weirs, error) { +func getWeirData(rm *RasModel, fn string, i int) (weirs, error) { weir := weirs{} f, err := rm.FileStore.GetObject(fn) @@ -503,7 +503,7 @@ func getWeirData(rm *RasModel, fn string, i int, shapeMap map[int]string) (weirs weir.NumGates++ case strings.HasPrefix(line, "IW Culv="): - conduit, err := getConduits(line, false, shapeMap) + conduit, err := getConduits(line, false) if err != nil { return weir, err } @@ -556,7 +556,7 @@ func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStruc case 2: var culvert culverts - culvert, i, err = getCulvertData(hsSc, i, data, conduitShapes) + culvert, i, err = getCulvertData(hsSc, i, data) if err != nil { return structures, err } @@ -573,7 +573,7 @@ func getHydraulicStructureData(rm *RasModel, fn string, idx int) (hydraulicStruc bData.NumBridges++ case 5: - weir, err := getWeirData(rm, fn, i, conduitShapes) + weir, err := getWeirData(rm, fn, i) if err != nil { return structures, err }