Working NEXRAD decode.

pull/58/head
Christopher Young 2015-09-22 19:40:29 -04:00
rodzic 0db905904f
commit f1ea6bd1a6
2 zmienionych plików z 229 dodań i 12 usunięć

Wyświetl plik

@ -0,0 +1,213 @@
package main
import (
"fmt"
"./uatparse"
"strconv"
"os"
"bufio"
)
// Most adapted from extract_nexrad.c
const (
BLOCK_WIDTH = float64(48.0/60.0)
WIDE_BLOCK_WIDTH = float64(96.0/60.0)
BLOCK_HEIGHT = float64(4.0/60.0)
BLOCK_THRESHOLD = 405000
BLOCKS_PER_RING = 450
)
type NEXRADFrame struct {
radar_type uint32
ts string
scale int
latNorth float64
lonWest float64
height float64
width float64
intensity []uint8 // Really only 4-bit values.
}
func block_location(block_num int, ns_flag bool, scale_factor int) (float64, float64, float64, float64) {
var realScale float64
if scale_factor == 1 {
realScale = float64(5.0)
} else if scale_factor == 2 {
realScale = float64(9.0)
} else {
realScale = float64(1.0)
}
if block_num >= BLOCK_THRESHOLD {
block_num = block_num & ^1
}
raw_lat := float64(BLOCK_HEIGHT * float64(int(float64(block_num) / float64(BLOCKS_PER_RING))))
raw_lon := float64(block_num % BLOCKS_PER_RING) * BLOCK_WIDTH
var lonSize float64
if block_num >= BLOCK_THRESHOLD {
lonSize = WIDE_BLOCK_WIDTH * realScale
} else {
lonSize = BLOCK_WIDTH * realScale
}
latSize := BLOCK_HEIGHT * realScale
if ns_flag { // Southern hemisphere.
raw_lat = 0 - raw_lat
} else {
raw_lat = raw_lat + BLOCK_HEIGHT
}
/*
if raw_lon > 180.0 {
raw_lon = raw_lon - 360.0
}*/
return raw_lat, raw_lon, latSize, lonSize
}
func decode_nexrad(f *uatparse.UATFrame) []NEXRADFrame {
ret := make([]NEXRADFrame, 0)
rle_flag := (uint32(f.FISB_data[0]) & 0x80) != 0
ns_flag := (uint32(f.FISB_data[0]) & 0x40) != 0
block_num := ((int(f.FISB_data[0]) & 0x0f) << 16) | (int(f.FISB_data[1]) << 8) | (int(f.FISB_data[2]))
scale_factor := (int(f.FISB_data[0]) & 0x30) >> 4
if rle_flag { // Single bin, RLE encoded.
lat, lon, h, w := block_location(block_num, ns_flag, scale_factor)
var tmp NEXRADFrame
tmp.radar_type = f.Product_id
tmp.ts = strconv.Itoa(int(f.FISB_hours)) + ":" + strconv.Itoa(int(f.FISB_minutes))
tmp.scale = scale_factor
tmp.latNorth = lat
tmp.lonWest = lon
tmp.height = h
tmp.width = w
tmp.intensity = make([]uint8, 0)
intensityData := f.FISB_data[3:]
for _, v := range intensityData {
intensity := uint8(v) & 0x7;
runlength := (uint8(v) >> 3) + 1
for runlength > 0 {
tmp.intensity = append(tmp.intensity, intensity)
runlength--
}
}
ret = append(ret, tmp)
} else {
var row_start int
var row_size int
if block_num >= 405000 {
row_start = block_num - ((block_num - 405000) % 225)
row_size = 225
} else {
row_start = block_num - (block_num % 450)
row_size = 450
}
row_offset := block_num - row_start
L := int(f.FISB_data[3] & 15)
for i := 0; i < L; i++ {
var bb int
if i == 0 {
bb = (int(f.FISB_data[3]) & 0xF0) | 0x08
} else {
bb = int(f.FISB_data[i+3])
}
for j := 0; j < 8; j++ {
if bb & (1 << uint(j)) != 0 {
row_x := (row_offset + 8*i + j - 3) % row_size
bn := row_start + row_x
lat, lon, h, w := block_location(bn, ns_flag, scale_factor)
var tmp NEXRADFrame
tmp.radar_type = f.Product_id
tmp.ts = strconv.Itoa(int(f.FISB_hours)) + ":" + strconv.Itoa(int(f.FISB_minutes))
tmp.scale = scale_factor
tmp.latNorth = lat
tmp.lonWest = lon
tmp.height = h
tmp.width = w
tmp.intensity = make([]uint8, 0)
for k := 0; k < 128; k++ {
z := uint8(0)
if f.Product_id == 64 {
z = 1
}
tmp.intensity = append(tmp.intensity, z)
}
ret = append(ret, tmp)
}
}
}
}
return ret
}
func parseInput(buf string) []NEXRADFrame {
ret := make([]NEXRADFrame, 0)
uatmsg, err := uatparse.New(buf)
if err != nil {
return ret
}
uatmsg.DecodeUplink()
for _, uatframe := range uatmsg.Frames {
if uatframe.Product_id != 63 && uatframe.Product_id != 64 { // It's neither CONUS nor Regional NEXRAD.
continue
}
p := decode_nexrad(uatframe)
if p != nil {
ret = append(ret, p...)
}
}
return ret
}
func main() {
if len(os.Args) < 2 {
fmt.Printf("%s <uat log>\n", os.Args[0])
return
}
fd, err := os.Open(os.Args[1])
if err != nil {
fmt.Printf("can't open file: %s\n", err.Error())
return
}
reader := bufio.NewReader(fd)
for {
buf, err := reader.ReadString('\n')
if err != nil {
fmt.Printf("lost stdin.\n")
break
}
z := parseInput(buf)
for _, zz := range z {
n := "Regional"
if zz.radar_type == 64 {
n = "CONUS"
}
fmt.Printf("NEXRAD %s %s %d %.0f %.0f %.0f %.0f ", n, zz.ts, zz.scale, zz.latNorth * 60, zz.lonWest * 60, zz.height * 60, zz.width * 60)
for _, intens := range zz.intensity {
fmt.Printf("%d", intens)
}
fmt.Printf("\n")
}
}
}

Wyświetl plik

@ -27,7 +27,10 @@ const (
)
type UATFrame struct {
raw_data []byte
Raw_data []byte
FISB_data []byte
FISB_hours uint32
FISB_minutes uint32
frame_length uint32
Frame_type uint32
@ -81,11 +84,7 @@ func dlac_decode(data []byte, data_len uint32) string {
func (f *UATFrame) decodeInfoFrame() {
f.Product_id = ((uint32(f.raw_data[0]) & 0x1f) << 6) | (uint32(f.raw_data[1]) >> 2)
if f.Product_id != 413 { // FIXME.
return
}
f.Product_id = ((uint32(f.Raw_data[0]) & 0x1f) << 6) | (uint32(f.Raw_data[1]) >> 2)
if f.Frame_type != 0 {
return // Not FIS-B.
@ -94,7 +93,7 @@ func (f *UATFrame) decodeInfoFrame() {
return // Too short for FIS-B.
}
t_opt := ((uint32(f.raw_data[1]) & 0x01) << 1) | (uint32(f.raw_data[2]) >> 7)
t_opt := ((uint32(f.Raw_data[1]) & 0x01) << 1) | (uint32(f.Raw_data[2]) >> 7)
if t_opt != 0 { //FIXME.
// fmt.Printf("don't know time format %d\n", t_opt)
@ -102,11 +101,16 @@ func (f *UATFrame) decodeInfoFrame() {
return
}
/* fisb_hours := (uint32(data[2]) & 0x7c) >> 2
fisb_minutes := ((uint32(data[2]) & 0x03) << 4) | (uint32(data[3]) >> 4)
*/
f.FISB_hours = (uint32(f.Raw_data[2]) & 0x7c) >> 2
f.FISB_minutes = ((uint32(f.Raw_data[2]) & 0x03) << 4) | (uint32(f.Raw_data[3]) >> 4)
fisb_length := f.frame_length - 4
fisb_data := f.raw_data[4:]
fisb_data := f.Raw_data[4:]
f.FISB_data = fisb_data
if f.Product_id != 413 { // Doesn't have text data that we'd be interested in, so we're done.
return
}
p := dlac_decode(fisb_data, fisb_length)
ret := make([]string, 0)
@ -178,7 +182,7 @@ func (u *UATMsg) DecodeUplink() error {
data = data[2 : frame_length+2]
thisFrame := new(UATFrame)
thisFrame.raw_data = data
thisFrame.Raw_data = data
thisFrame.frame_length = frame_length
thisFrame.Frame_type = frame_type