Guppi RAW output

Makefile

@@ -39,6 +39,39 @@ LINKH5:= $(LIBDIR_H5) -l $(LIBHDF5) -l $(LIBHDF5_HL)
 
 # End HDF5 definitions
 
+# Begin GuppiRawC99 definitions
+
+# if we have pkg-config installed, do we have guppirawc99 libraries installed?
+ifeq ($(HAVEPKG),yes)
+  HAVEGRC99 = $(shell pkg-config --exists guppirawc99 && echo yes)
+endif
+
+# depending on the above questions use, in order of priority:
+# 1. User supplied values for INCDIR_GRC99 and LIBDIR_GRC99
+# 2. Whatever pkg-config supplies for those variables 
+# 3. Failing either the above, default values /usr/local/include and /usr/local/lib
+
+ifeq ($(INCDIR_GRC99),)
+  ifeq ($(HAVEGRC99),yes)
+    INCDIR_GRC99 = $(shell pkg-config --cflags guppirawc99)
+  else
+    INCDIR_GRC99 = -I/usr/local/include
+  endif
+endif
+
+ifeq ($(LIBDIR_GRC99),)
+  ifeq ($(HAVEGRC99),yes)
+    LIBDIR_GRC99 = $(shell pkg-config --libs guppirawc99)
+  else
+    LIBDIR_GRC99 = -L/usr/local/lib
+  endif
+endif
+
+LIBGRC99= :libguppirawc99.so
+LINKGRC99:= $(LIBDIR_GRC99) -l $(LIBGRC99)
+
+# End GuppiRawC99 definitions
+
 CUDA_DIR ?= $(CUDA_ROOT)
 CUDA_PATH ?= $(CUDA_DIR)
 
@@ -47,7 +80,7 @@ CXX           ?= g++
 HOST_COMPILER ?= $(CXX)
 NVCC          := $(CUDA_PATH)/bin/nvcc -ccbin $(HOST_COMPILER)
 
-CFLAGS = -ggdb -fPIC -I$(CUDA_PATH)/include $(INCDIR_H5)
+CFLAGS = -ggdb -fPIC -I$(CUDA_PATH)/include $(INCDIR_H5) $(INCDIR_GRC99)
 ifdef DEBUG_CALLBACKS
 CFLAGS += -DDEBUG_CALLBACKS=$(DEBUG_CALLBACKS)
 endif
@@ -98,7 +131,7 @@ rawspec_gpu.o: rawspec.h rawspec_version.h cufft_error_name.h
 rawspec_socket.o: rawspec_socket.h rawspec.h \
                   rawspec_callback.h rawspec_fbutils.h
 rawspectest.o: rawspec.h
-rawspec_rawutils.o: rawspec_rawutils.h hget.h
+rawspec_rawutils.o: rawspec_rawutils.h
 
 # Begin fbh5 objects
 fbh5_open.o: fbh5_defs.h rawspec_callback.h rawspec_fbutils.h
@@ -111,11 +144,11 @@ fbh5_util.o: fbh5_defs.h rawspec_callback.h rawspec_fbutils.h
 	$(VERBOSE) $(NVCC) $(NVCC_FLAGS) -dc $(GENCODE_FLAGS) -o $@ -c $<
 
 librawspec.so: rawspec_gpu.o rawspec_fbutils.o rawspec_rawutils.o fbh5_open.o fbh5_close.o fbh5_write.o fbh5_util.o
-	$(VERBOSE) $(NVCC) -shared $(NVCC_FLAGS) $(GENCODE_FLAGS) -o $@ $^ $(CUDA_STATIC_LIBS) $(LINKH5)
+	$(VERBOSE) $(NVCC) -shared $(NVCC_FLAGS) $(GENCODE_FLAGS) -o $@ $^ $(CUDA_STATIC_LIBS) $(LINKH5) $(LINKGRC99)
 
 rawspec: librawspec.so
 rawspec: rawspec.o rawspec_file.o rawspec_socket.o 
-	$(VERBOSE) $(NVCC) $(NVCC_FLAGS) $(GENCODE_FLAGS) -o $@ $^ -L. -lrawspec $(LINKH5)
+	$(VERBOSE) $(NVCC) $(NVCC_FLAGS) $(GENCODE_FLAGS) -o $@ $^ -L. -lrawspec $(LINKH5) $(LINKGRC99)
 
 rawspectest: librawspec.so
 rawspectest: rawspectest.o

compare.jl

@@ -0,0 +1,77 @@
+using Printf
+using FFTW
+using Blio: GuppiRaw
+
+struct Result
+  message::Union{String, Nothing}
+  value::Bool
+end
+
+function _fft(guppidata, points)
+  data = reshape(guppidata, (size(guppidata, 1), points, :, size(guppidata)[3:end]...))
+  # [pol, fine_spec, time, coarse_chan, antenna]
+  spectra = fft(data, 2)
+  # [pol, fine_chan, time, antenna]
+  spectra = cat((spectra[:, :, :, i, :] for i in 1:size(spectra,4))..., dims=2)
+  # [pol, time, fine_chan, antenna]
+  return permutedims(spectra, [1,3,2,4])
+end
+
+function accumulate(guppidata, samples)
+  data = reshape(guppidata, (size(guppidata, 1), samples, :, size(guppidata)[3:end]...))
+  acc = sum(data, dims=2)
+  return reshape(acc, (size(guppidata, 1), :, size(guppidata)[3:end]...))
+end
+
+function mapToFloat(value::Integer, type::Type)
+  return value < 0 ? -1.0*(value / typemin(type)) : value / typemax(type)
+end
+
+function mapToFloat(value::Complex{<:Integer}, type::Type)
+  return complex(mapToFloat(real(value), real(type)), mapToFloat(imag(value), real(type)))
+end
+
+function compare(i_data, o_data, atol=0.01)::Result
+  if size(i_data) != size(o_data)
+    return Result(@sprintf("Shape mismatch: %s != %s", size(i_data), size(o_data)), false)
+  end
+  dims_correct = Array{Bool}(undef, size(i_data)[2:end])
+  for i in CartesianIndices(dims_correct)
+    dims_correct[i] = all(isapprox.(real(i_data[:, i]), real(o_data[:, i]), atol=atol)) && all(isapprox.(imag(i_data[:, i]), imag(o_data[:, i]), atol=atol))
+    if !dims_correct[i]
+      println(Result(@sprintf("Pol data mismatch @ %s: %s != %s\n\t(atol: %s)", i, i_data[:, i], o_data[:, i], i_data[:, i] - o_data[:, i]), false))
+    end
+  end
+
+  Result(nothing, all(dims_correct))
+end
+
+i_grheader = GuppiRaw.Header()
+o_grheader = GuppiRaw.Header()
+
+i_fio = open("/mnt/buf0/test.0000.raw", "r")
+o_fio = open("/mnt/buf0/test.rawspec.0000.raw", "r")
+
+  read!(i_fio, i_grheader)
+  i_data = Array(i_grheader)
+  read!(i_fio, i_data)
+  i_data = mapToFloat.(i_data, eltype(i_data))
+
+  read!(o_fio, o_grheader)
+  o_data = Array(o_grheader)
+  read!(o_fio, o_data)
+
+  fftpoints = div(size(o_data, 3), size(i_data, 3))
+  accumulation = div(div(size(i_data, 2), size(o_data, 2)), fftpoints)
+  println("fftpoints: ", fftpoints)
+  println("accumulation: ", accumulation)
+
+  upchannelized = fftpoints == 1 ? i_data : _fft(i_data, fftpoints)
+  accum_fine = accumulation == 1 ? upchannelized : accumulate(upchannelized, accumulation)
+
+  atol = 0.015 * max(log2(fftpoints), fftpoints) * accumulation
+
+  println(compare(accum_fine, o_data, atol))
+
+close(i_fio)
+close(o_fio)