forked from jda/srtm
-
Notifications
You must be signed in to change notification settings - Fork 0
/
file.go
129 lines (111 loc) · 3.34 KB
/
file.go
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
package srtm
import (
"compress/gzip"
"encoding/binary"
"fmt"
"io/ioutil"
"os"
"regexp"
"strings"
"github.com/pkg/errors"
)
// ErrInvalidHGTFileName is returned when a HGT file name does not match
// pattern of a valid file which indicates lat/lon
var ErrInvalidHGTFileName = errors.New("invalid HGT file name")
var srtmParseName = regexp.MustCompile(`(N|S)(\d\d)(E|W)(\d\d\d)\.hgt(\.gz)?`)
// ReadFile is a helper func around Read that reads a SRTM file, decompressing
// if necessary, and returns SRTM elevation data
func ReadFile(file string) (sw *LatLng, squareSize int, elevations []int16, err error) {
f, err := os.Open(file)
if err != nil {
return sw, squareSize, elevations, err
}
defer f.Close()
if strings.HasSuffix(file, ".gz") {
rdr, err := gzip.NewReader(f)
if err != nil {
return sw, squareSize, elevations, err
}
defer rdr.Close()
bytes, err := ioutil.ReadAll(rdr)
if err != nil {
return sw, squareSize, elevations, err
}
return Read(file, bytes)
}
bytes, err := ioutil.ReadAll(f)
if err != nil {
return sw, squareSize, elevations, err
}
return Read(file, bytes)
}
// Read reads elevation for points from a SRTM file
func Read(fname string, bytes []byte) (sw *LatLng, squareSize int, elevations []int16, err error) {
if len(bytes) == 12967201*2 {
// 1 arcsecond
squareSize = 3601
} else if len(bytes) == 1442401*2 {
// 3 arcseconds
squareSize = 1201
} else {
return sw, squareSize, elevations, fmt.Errorf("hgt file cannot identified (only 1 arcsecond and 3 arcsecond supported, file size = %d)", len(bytes))
}
sw, err = southWest(fname)
if err != nil {
return sw, squareSize, elevations, errors.Wrap(err, "could not get corner coordinates from file name")
}
elevations = make([]int16, squareSize*squareSize)
// Latitude
for row := 0; row < squareSize; row++ {
// Longitude
for col := 0; col < squareSize; col++ {
idx := row*squareSize + col
elevations[idx] = int16(binary.BigEndian.Uint16(bytes[idx*2 : idx*2+2]))
}
}
return sw, squareSize, elevations, nil
}
// Meta reads meta information of SRTM file
func Meta(fname string, size int64) (sw *LatLng, squareSize int, err error) {
if size == 12967201*2 {
// 1 arcsecond
squareSize = 3601
} else if size == 1442401*2 {
// 3 arcseconds
squareSize = 1201
} else {
return sw, squareSize, fmt.Errorf("hgt file cannot identified (only 1 arcsecond and 3 arcsecond supported, file size = %d)", size)
}
sw, err = southWest(fname)
return sw, squareSize, err
}
// sw returns the southwest point contained in a HGT file.
// Coordinates in the file are relative to this point
func southWest(file string) (p *LatLng, err error) {
fnameParts := srtmParseName.FindStringSubmatch(file)
if fnameParts == nil {
return p, ErrInvalidHGTFileName
}
swLatitude, err := dToDecimal(fnameParts[1] + fnameParts[2])
if err != nil {
return p, errors.Wrap(err, "could not get Latitude from file name")
}
swLongitude, err := dToDecimal(fnameParts[3] + fnameParts[4])
if err != nil {
return p, errors.Wrap(err, "could not get Longitude from file name")
}
return &LatLng{
Latitude: swLatitude,
Longitude: swLongitude,
}, err
}
// IsHGT returns true if fname appears to be a SRTM HGT file
func IsHGT(fname string) bool {
if strings.HasSuffix(fname, ".hgt") {
return true
}
if strings.HasSuffix(fname, ".hgt.gz") {
return true
}
return false
}