packages feed

sdr 0.1.0.5 → 0.1.0.6

raw patch · 14 files changed

+445/−80 lines, 14 filesdep +mwc-randomdep ~Chartdep ~Chart-cairodep ~cerealPVP: major bump suggested

API removals or changes: PVP suggests a major version bump

Dependencies added: mwc-random

Dependency ranges changed: Chart, Chart-cairo, cereal, optparse-applicative, sdr, vector

API changes (from Hackage documentation)

- SDR.FFT: fftFixup :: (Vector v n, Num n) => Int -> v n
+ SDR.Filter: fastResamplerAVXC :: Int -> Int -> [Float] -> IO (Resampler IO Vector MVector (Complex Float))
+ SDR.Filter: fastResamplerC :: CPUInfo -> Int -> Int -> [Float] -> IO (Resampler IO Vector MVector (Complex Float))
+ SDR.Filter: fastResamplerCC :: Int -> Int -> [Float] -> IO (Resampler IO Vector MVector (Complex Float))
+ SDR.Filter: fastResamplerSSEC :: Int -> Int -> [Float] -> IO (Resampler IO Vector MVector (Complex Float))
+ SDR.FilterInternal: mkResamplerC :: ResampleR -> Int -> Int -> Int -> [Float] -> IO (Int -> Int -> Vector (Complex Float) -> MVector RealWorld (Complex Float) -> IO Int)
+ SDR.FilterInternal: resample2RC_c :: CInt -> CInt -> CInt -> CInt -> Ptr CInt -> Ptr (Ptr CFloat) -> Ptr CFloat -> Ptr CFloat -> IO CInt
+ SDR.FilterInternal: resampleAVXRC_c :: CInt -> CInt -> CInt -> CInt -> Ptr CInt -> Ptr (Ptr CFloat) -> Ptr CFloat -> Ptr CFloat -> IO CInt
+ SDR.FilterInternal: resampleCAVXRC :: ResampleRC
+ SDR.FilterInternal: resampleCRC :: ResampleRC
+ SDR.FilterInternal: resampleCSSERC :: ResampleRC
+ SDR.FilterInternal: resampleCSSERC_c :: CInt -> CInt -> CInt -> CInt -> Ptr CInt -> Ptr (Ptr CFloat) -> Ptr CFloat -> Ptr CFloat -> IO CInt
+ SDR.FilterInternal: resampleFFIC :: (Ptr CFloat -> Ptr CFloat -> IO CInt) -> Vector (Complex Float) -> MVector RealWorld (Complex Float) -> IO Int
+ SDR.FilterInternal: type ResampleRC = Int -> Int -> [Float] -> IO (Int -> Int -> Vector (Complex Float) -> MVector RealWorld (Complex Float) -> IO Int)
+ SDR.Util: agc :: (Num a, Storable a, RealFloat a) => a -> a -> a -> Vector (Complex a) -> (a, Vector (Complex a))
+ SDR.Util: agcPipe :: (Num a, Storable a, RealFloat a, Monad m) => a -> a -> Pipe (Vector (Complex a)) (Vector (Complex a)) m ()
+ SDR.Util: complexFloatToInterleavedIQSigned2048 :: (Integral b, RealFrac a, Vector v1 (Complex a), Vector v2 b) => v1 (Complex a) -> v2 b
+ SDR.Util: complexFloatToInterleavedIQSignedWord :: Vector (Complex Float) -> Vector CShort
+ SDR.Util: halfBandUp :: (Vector v n, Num n) => Int -> v n
+ SDR.Util: interleavedIQSignedWordToFloatFast :: CPUInfo -> Vector CShort -> Vector (Complex Float)
+ SDR.Util: streamRandom :: PrimMonad m => Int -> Producer (Vector Float) m ()
+ SDR.Util: streamString :: (FiniteBits b, Monad m) => [b] -> Int -> Producer (Vector Float) m ()
+ SDR.VectorUtils: copyInto :: (PrimMonad m, MVector vm a, Vector v a) => vm (PrimState m) a -> v a -> m ()
+ SDR.VectorUtils: vUnfoldrM :: (PrimMonad m, Vector v x) => Int -> (acc -> m (x, acc)) -> acc -> m (v x, acc)
- SDR.Demod: fmDemodStr :: (RealFloat a) => Complex a -> Stream (Complex a) -> Stream a
+ SDR.Demod: fmDemodStr :: (RealFloat a, Monad m) => Complex a -> Stream m (Complex a) -> Stream m a
- SDR.VectorUtils: fill :: (PrimMonad m, Functor m, MVector vm a) => MStream m a -> vm (PrimState m) a -> m ()
+ SDR.VectorUtils: fill :: (PrimMonad m, Functor m, MVector vm a) => Bundle v a -> vm (PrimState m) a -> m ()

Files

benchmarks/Benchmarks.hs view
@@ -71,6 +71,10 @@     resampler4 <- resampleCSSERR interpolation decimation coeffsList     resampler5 <- resampleCAVXRR interpolation decimation coeffsList +    resampler3C <- resampleCRC    interpolation decimation coeffsList+    resampler4C <- resampleCSSERC interpolation decimation coeffsList+    resampler5C <- resampleCAVXRC interpolation decimation coeffsList+     --Benchmarks     defaultMain [             bgroup "filter" [@@ -124,7 +128,10 @@                     (hasAVX,     bench "cAVX"        $ nfIO $ resampler5               (num `quot` decimation) 0 inBuf outBuf)                 ],                 bgroup "complex" $ hasFeatures [-                    (const True, bench "highLevel"   $ nfIO $ resampleHighLevel        interpolation decimation coeffs 0 (num `quot` decimation) inBufComplex outBufComplex)+                    (const True, bench "highLevel"   $ nfIO $ resampleHighLevel        interpolation decimation coeffs 0 (num `quot` decimation) inBufComplex outBufComplex),+                    (const True, bench "c"           $ nfIO $ resampler3C               (num `quot` decimation) 0 inBufComplex outBufComplex),+                    (const True, bench "SSE"         $ nfIO $ resampler4C               (num `quot` decimation) 0 inBufComplex outBufComplex),+                    (const True, bench "AVX"         $ nfIO $ resampler5C               (num `quot` decimation) 0 inBufComplex outBufComplex)                 ]             ],             bgroup "scaling" $ hasFeatures [
c_sources/convert.c view
@@ -84,3 +84,19 @@     } } +void convertBladeRFTransmit(int num, float *in, int16_t *out){+    int i;+    for(i=0; i<num; i++){+        float val   = in[i];        // [-1, 1)+        val         = val + 1;      // [0, 2)+        val         = val * 2048;   // [0, 4096)+        int16_t res = (int16_t)val; // [0, 4095]+        res         = res - 2048;   // [-2048, 2047]++        if(res > 2047)  res = 2047;+        if(res < -2048) res = -2048;++        out[i]      = res;+    }+}+
c_sources/resample.c view
@@ -13,7 +13,7 @@ /*  * Rational downsampling  */-void resample(int buf_size, int coeff_size, int interpolation, int decimation, int filter_offset, float *coeffs, float *in_buf, float *out_buf){+void resampleRR(int buf_size, int coeff_size, int interpolation, int decimation, int filter_offset, float *coeffs, float *in_buf, float *out_buf){     int j, k, l;     int input_offset = 0;     for(k=0; k<buf_size; k++) {@@ -31,7 +31,7 @@     } } -int resample2(int buf_size, int num_coeffs, int starting_group, int num_groups, int *increments, float **coeffs, float *in_buf, float *out_buf){+int resample2RR(int buf_size, int num_coeffs, int starting_group, int num_groups, int *increments, float **coeffs, float *in_buf, float *out_buf){     int   i, j;     int   group      = starting_group;     float *start_ptr = in_buf;@@ -78,6 +78,61 @@         _mm_store_ss(out_buf + i, res);          startPtr  += increments[group];+        group++;+        if(group == num_groups) +            group = 0;+    }++    return group;+}++int resample2RC(int buf_size, int num_coeffs, int starting_group, int num_groups, int *increments, float **coeffs, float *in_buf, float *out_buf){+    int   i;+    int   group      = starting_group;+    float *start_ptr = in_buf;++    for(i=0; i<buf_size*2; i+=2){+        dotprod_C(num_coeffs, coeffs[group], start_ptr, out_buf + i);++        start_ptr  += (2 * increments[group]);+        group++;+        if(group == num_groups) +            group = 0;+    }++    return group;+}++int resampleSSERC(int buf_size, int num_coeffs, int starting_group, int num_groups, int *increments, float **coeffs, float *in_buf, float *out_buf){+    int   i;+    int   group      = starting_group;+    float *startPtr  = in_buf;++    for(i=0; i<buf_size*2; i+=2){+        __m128 accum = sse_dotprod_C(num_coeffs, coeffs[group], startPtr);+        accum = sse_hadd_C(accum);+        store_complex(out_buf + i, accum);++        startPtr  += (2 * increments[group]);+        group++;+        if(group == num_groups) +            group = 0;+    }++    return group;+}++int resampleAVXRC(int buf_size, int num_coeffs, int starting_group, int num_groups, int *increments, float **coeffs, float *in_buf, float *out_buf){+    int   i;+    int   group      = starting_group;+    float *startPtr  = in_buf;++    for(i=0; i<buf_size*2; i+=2){+        __m256 accum = avx_dotprod_C(num_coeffs, coeffs[group], startPtr);+        __m128 res   = avx_hadd_C(accum);+        store_complex(out_buf + i, res);++        startPtr  += (2 * increments[group]);         group++;         if(group == num_groups)              group = 0;
hs_sources/SDR/Demod.hs view
@@ -8,18 +8,20 @@     ) where  import Data.Complex-import Data.Vector.Generic       as VG-import Data.Vector.Fusion.Stream+import qualified Data.Vector.Generic               as VG+import qualified Data.Vector.Fusion.Stream.Monadic as VFSM+import qualified Data.Vector.Fusion.Bundle.Monadic as VFBM+import qualified Data.Vector.Fusion.Bundle.Size    as VFBS  import SDR.VectorUtils import Pipes  -- | FM demodulate a stream of complex samples {-# INLINE fmDemodStr #-}-fmDemodStr :: (RealFloat a) +fmDemodStr :: (RealFloat a, Monad m)             => Complex a          -- ^ The starting sample - i.e. the last sample in the last buffer-           -> Stream (Complex a) -- ^ The input stream-           -> Stream a           -- ^ The output stream+           -> VFSM.Stream m (Complex a) -- ^ The input stream+           -> VFSM.Stream m a           -- ^ The output stream fmDemodStr = mapAccumMV func      where     {-# INLINE func #-}@@ -27,15 +29,15 @@  -- | FM demodulate a vector of complex samples {-# INLINE fmDemodVec #-}-fmDemodVec :: (RealFloat a, Vector v (Complex a), Vector v a) +fmDemodVec :: (RealFloat a, VG.Vector v (Complex a), VG.Vector v a)             => Complex a     -- ^ The starting sample - i.e. the last sample in the last buffer            -> v (Complex a) -- ^ The input Vector            -> v a           -- ^ The output Vector-fmDemodVec init = unstream . fmDemodStr init . stream+fmDemodVec init inp = VG.unstream $ flip VFBM.fromStream (VFBS.Exact $ VG.length inp) $ fmDemodStr init $ VFBM.elements $ VG.stream inp  -- | Pipe that performs FM demodulation {-# INLINE fmDemod #-}-fmDemod :: (RealFloat a, Vector v (Complex a), Vector v a) => Pipe (v (Complex a)) (v a) IO ()+fmDemod :: (RealFloat a, VG.Vector v (Complex a), VG.Vector v a) => Pipe (v (Complex a)) (v a) IO () fmDemod = func 0     where     func lastSample = do
hs_sources/SDR/FFT.hs view
@@ -3,7 +3,6 @@ {-| Fast FFTs using FFTW -} module SDR.FFT (     -- * Windows-    fftFixup,     hamming,     hanning,     blackman,@@ -23,32 +22,21 @@ import           Control.Concurrent           hiding (yield) import qualified Data.Map                     as Map -import           Data.Vector.Generic          as VG-import           Data.Vector.Generic.Mutable  as VGM-import           Data.Vector.Storable         as VS-import           Data.Vector.Storable.Mutable as VSM-import           Data.Vector.Fusion.Stream    as VFS+import qualified Data.Vector.Generic          as VG+import qualified Data.Vector.Storable         as VS+import qualified Data.Vector.Storable.Mutable as VSM  import           Pipes import           Numeric.FFTW  import           SDR.FilterDesign+import           SDR.VectorUtils  mallocForeignBufferAligned :: forall a. Storable a => Int -> IO (ForeignPtr a) mallocForeignBufferAligned elems = do     ptr <- fftwMalloc $ fromIntegral $ elems * sizeOf (undefined :: a)     newForeignPtr fftwFreePtr ptr --- | Compute a vector of alternating 1s and 0s of the given size.-fftFixup :: (VG.Vector v n, Num n) -         => Int -- ^ The length of the Vector-         -> v n-fftFixup size = VG.generate size func-    where-    func idx -        | even idx  = 1-        | otherwise = -1- -- | Creates a Pipe that performs a complex to complex DFT. fftw :: (VG.Vector v (Complex CDouble))       => Int -- ^ The size of the input and output buffers@@ -66,7 +54,7 @@         ina <- lift $ mallocForeignBufferAligned samples         let inv = VSM.unsafeFromForeignPtr0 ina samples -        lift $ VGM.fill inv $ VFS.mapM return $ VG.stream inv'+        lift $ copyInto inv inv'          let (fp, offset, length) = VSM.unsafeToForeignPtr inv @@ -94,7 +82,7 @@         ina <- lift $ mallocForeignBufferAligned samples         let inv = VSM.unsafeFromForeignPtr0 ina samples -        lift $ VGM.fill inv $ VFS.mapM return $ VG.stream inv'+        lift $ copyInto inv inv'         let (fp, offset, length) = VSM.unsafeToForeignPtr inv          lift $ withForeignPtr fp  $ \fpp -> @@ -132,7 +120,7 @@         ina <- mallocForeignBufferAligned samples         let inv = VSM.unsafeFromForeignPtr0 ina samples -        VGM.fill inv $ VFS.mapM return $ VG.stream res+        copyInto inv res          let (fp, offset, length) = VSM.unsafeToForeignPtr inv 
hs_sources/SDR/Filter.hs view
@@ -76,6 +76,12 @@     fastResamplerAVXR,     fastResamplerR, +    -- *** Complex Data+    fastResamplerCC,+    fastResamplerSSEC,+    fastResamplerAVXC,+    fastResamplerC,+     -- * Filter     firFilter, @@ -100,7 +106,7 @@ import           Pipes  import           SDR.Util-import           SDR.FilterInternal          hiding (mkResampler)+import           SDR.FilterInternal          hiding (mkResampler, mkResamplerC) import           SDR.CPUID  {- | A `Filter` contains all of the information needed by the `filterr` @@ -418,6 +424,25 @@         startDat                = (0, 0)     return $ Resampler {..} +mkResamplerC :: Int+             -> ResampleRC+             -> Int+             -> Int+             -> [Float] +             -> IO (Resampler IO VS.Vector VS.MVector (Complex Float)) +mkResamplerC sizeMultiple filterFunc interpolationR decimationR coeffs = do+    let vCoeffs     = VG.fromList coeffs+    evaluate vCoeffs+    resamp <- filterFunc interpolationR decimationR coeffs+    let resampleOne   v w x y   = func1 <$> resamp (fst v) w x y+        resampleCross (group, offset) count x y z = do +            offset' <- resampleCrossHighLevel interpolationR decimationR vCoeffs offset count x y z+            return (((group + count) `mod` interpolationR, offset'), offset')+        numCoeffsR              = roundUp (length coeffs) (interpolationR * sizeMultiple)+        func1 group             = let offset = interpolationR - 1 - ((interpolationR + group * decimationR - 1) `mod` interpolationR) in ((group, offset), offset)+        startDat                = (0, 0)+    return $ Resampler {..}+ -- | Returns a fast Resampler data structure implemented in C. For filtering real data with real coefficients. fastResamplerCR :: Int                                          -- ^ The interpolation factor                 -> Int                                          -- ^ The decimation factor@@ -446,6 +471,35 @@                -> [Float]                                      -- ^ The filter coefficients                -> IO (Resampler IO VS.Vector VS.MVector Float) -- ^ The `Resampler` data structure fastResamplerR info = featureSelect info fastResamplerCR [(hasAVX, fastResamplerAVXR), (hasSSE42, fastResamplerSSER)]++-- | Returns a fast Resampler data structure implemented in C. For filtering complex data with real coefficients.+fastResamplerCC :: Int                                          -- ^ The interpolation factor+                -> Int                                          -- ^ The decimation factor+                -> [Float]                                      -- ^ The filter coefficients+                -> IO (Resampler IO VS.Vector VS.MVector (Complex Float)) -- ^ The `Resampler` data structure+fastResamplerCC = mkResamplerC 1 resampleCRC++-- | Returns a fast Resampler data structure implemented in C using SSE instructions. For filtering complex data with real coefficients.+fastResamplerSSEC :: Int                                          -- ^ The interpolation factor+                  -> Int                                          -- ^ The decimation factor+                  -> [Float]                                      -- ^ The filter coefficients+                  -> IO (Resampler IO VS.Vector VS.MVector (Complex Float)) -- ^ The `Resampler` data structure+fastResamplerSSEC = mkResamplerC 4 resampleCSSERC++-- | Returns a fast Resampler data structure implemented in C using AVX instructions. For filtering complex data with real coefficients.+fastResamplerAVXC :: Int                                          -- ^ The interpolation factor+                  -> Int                                          -- ^ The decimation factor+                  -> [Float]                                      -- ^ The filter coefficients+                  -> IO (Resampler IO VS.Vector VS.MVector (Complex Float)) -- ^ The `Resampler` data structure+fastResamplerAVXC = mkResamplerC 8 resampleCAVXRC++-- | Returns a fast Resampler data structure implemented in C using the fastest SIMD instruction set your processor supports. For resampling complex data with real coefficients.+fastResamplerC :: CPUInfo                                      -- ^ The CPU's capabilities+               -> Int                                          -- ^ The interpolation factor+               -> Int                                          -- ^ The decimation factor+               -> [Float]                                      -- ^ The filter coefficients+               -> IO (Resampler IO VS.Vector VS.MVector (Complex Float)) -- ^ The `Resampler` data structure+fastResamplerC info = featureSelect info fastResamplerCC [(hasAVX, fastResamplerAVXC), (hasSSE42, fastResamplerSSEC)]  data Buffer v a = Buffer {     buffer :: v a,
hs_sources/SDR/FilterDesign.hs view
@@ -80,7 +80,11 @@  --ts is really the ratio ts / sampling period -- | Square root raised cosine-srrc :: (Ord a, Floating a) => Int -> Int -> a -> [a]+srrc :: (Ord a, Floating a) +     => Int -- ^ size: from [-n .. n]+     -> Int -- ^ sampling period+     -> a   -- ^ beta+     -> [a] srrc n ts beta = map func [(-n) .. n]     where     func x 
hs_sources/SDR/FilterInternal.hs view
@@ -20,15 +20,14 @@ import qualified Data.Vector.Generic.Mutable       as VGM import qualified Data.Vector.Storable              as VS import qualified Data.Vector.Storable.Mutable      as VSM-import qualified Data.Vector.Fusion.Stream         as VFS-import qualified Data.Vector.Fusion.Stream.Monadic as VFSM+import qualified Data.Vector.Fusion.Bundle         as VFB  import           SDR.VectorUtils import           SDR.Util  {-# INLINE filterHighLevel #-} filterHighLevel :: (PrimMonad m, Functor m, Num a, Mult a b, VG.Vector v a, VG.Vector v b, VGM.MVector vm a) => v b -> Int -> v a -> vm (PrimState m) a -> m ()-filterHighLevel coeffs num inBuf outBuf = fill (VFSM.generate num dotProd) outBuf+filterHighLevel coeffs num inBuf outBuf = fill (VFB.generate num dotProd) outBuf     where     dotProd offset = VG.sum $ VG.zipWith mult (VG.unsafeDrop offset inBuf) coeffs @@ -156,7 +155,7 @@ decimateHighLevel :: (PrimMonad m, Functor m, Num a, Mult a b, VG.Vector v a, VG.Vector v b, VGM.MVector vm a) => Int -> v b -> Int -> v a -> vm (PrimState m) a -> m () decimateHighLevel factor coeffs num inBuf outBuf = fill x outBuf     where -    x = VFSM.map dotProd (VFSM.iterateN num (+ factor) 0)+    x = VFB.map dotProd (VFB.iterateN num (+ factor) 0)     dotProd offset = VG.sum $ VG.zipWith mult (VG.unsafeDrop offset inBuf) coeffs  type DecimateCRR = CInt -> CInt -> CInt -> Ptr CFloat -> Ptr CFloat -> Ptr CFloat -> IO ()@@ -265,7 +264,7 @@         | otherwise = return filterOffset     dotProd filterOffset offset = VG.sum $ VG.zipWith mult (VG.unsafeDrop offset inBuf) (stride interpolation (VG.unsafeDrop filterOffset coeffs)) -foreign import ccall unsafe "resample"+foreign import ccall unsafe "resampleRR"     resample_c :: CInt -> CInt -> CInt -> CInt -> CInt -> Ptr CFloat -> Ptr CFloat -> Ptr CFloat -> IO ()  resampleCRR :: Int -> Int -> Int -> Int -> VS.Vector Float -> VS.Vector Float -> VS.MVector RealWorld Float -> IO ()@@ -325,6 +324,12 @@         VS.unsafeWith (unsafeCoerce outBuf) $ \oPtr ->              func iPtr oPtr +resampleFFIC :: (Ptr CFloat -> Ptr CFloat -> IO CInt) -> VS.Vector (Complex Float) -> VSM.MVector RealWorld (Complex Float) -> IO Int+resampleFFIC func inBuf outBuf = liftM fromIntegral $+    VS.unsafeWith (unsafeCoerce inBuf) $ \iPtr -> +        VS.unsafeWith (unsafeCoerce outBuf) $ \oPtr -> +            func iPtr oPtr+ type ResampleR = CInt -> CInt -> CInt -> CInt -> Ptr CInt -> Ptr (Ptr CFloat) -> Ptr CFloat -> Ptr CFloat -> IO CInt  mkResampler :: ResampleR -> Int -> Int -> Int -> [Float] -> IO (Int -> Int -> VS.Vector Float -> VS.MVector RealWorld Float -> IO Int)@@ -338,7 +343,7 @@  type ResampleRR = Int -> Int -> [Float] -> IO (Int -> Int -> VS.Vector Float -> VS.MVector RealWorld Float -> IO Int) -foreign import ccall unsafe "resample2"+foreign import ccall unsafe "resample2RR"     resample2_c :: CInt -> CInt -> CInt -> CInt -> Ptr CInt -> Ptr (Ptr CFloat) -> Ptr CFloat -> Ptr CFloat -> IO CInt  resampleCRR2 :: ResampleRR@@ -356,6 +361,35 @@ resampleCAVXRR :: ResampleRR resampleCAVXRR = mkResampler resampleAVXRR_c 8 +type ResampleRC = Int -> Int -> [Float] -> IO (Int -> Int -> VS.Vector (Complex Float) -> VS.MVector RealWorld (Complex Float) -> IO Int)++mkResamplerC :: ResampleR -> Int -> Int -> Int -> [Float] -> IO (Int -> Int -> VS.Vector (Complex Float) -> VS.MVector RealWorld (Complex Float) -> IO Int)+mkResamplerC func n interpolation decimation coeffs = do+    groupsP     <- mapM newArray $ map (map realToFrac) groups+    groupsPP    <- newArray groupsP+    incrementsP <- newArray $ map fromIntegral increments+    return $ \offset num -> resampleFFIC $ func (fromIntegral num) (fromIntegral numCoeffs) (fromIntegral offset) (fromIntegral numGroups) incrementsP groupsPP+    where+    Coeffs {..} = prepareCoeffs n interpolation decimation coeffs++foreign import ccall unsafe "resample2RC"+    resample2RC_c :: CInt -> CInt -> CInt -> CInt -> Ptr CInt -> Ptr (Ptr CFloat) -> Ptr CFloat -> Ptr CFloat -> IO CInt++resampleCRC :: ResampleRC+resampleCRC = mkResamplerC resample2RC_c 1++foreign import ccall unsafe "resampleSSERC"+    resampleCSSERC_c :: CInt -> CInt -> CInt -> CInt -> Ptr CInt -> Ptr (Ptr CFloat) -> Ptr CFloat -> Ptr CFloat -> IO CInt++resampleCSSERC :: ResampleRC+resampleCSSERC = mkResamplerC resampleCSSERC_c 4++foreign import ccall unsafe "resampleAVXRC"+    resampleAVXRC_c :: CInt -> CInt -> CInt -> CInt -> Ptr CInt -> Ptr (Ptr CFloat) -> Ptr CFloat -> Ptr CFloat -> IO CInt++resampleCAVXRC :: ResampleRC+resampleCAVXRC = mkResamplerC resampleAVXRC_c 8+ {-  - Cross buffer -}@@ -364,12 +398,12 @@ decimateCrossHighLevel :: (PrimMonad m, Functor m, Num a, Mult a b, VG.Vector v a, VG.Vector v b, VGM.MVector vm a) => Int -> v b -> Int -> v a -> v a -> vm (PrimState m) a -> m () decimateCrossHighLevel factor coeffs num lastBuf nextBuf outBuf = fill x outBuf     where-    x = VFSM.map dotProd (VFSM.iterateN num (+ factor) 0)+    x = VFB.map dotProd (VFB.iterateN num (+ factor) 0)     dotProd i = VG.sum $ VG.zipWith mult (VG.unsafeDrop i lastBuf VG.++ nextBuf) coeffs  {-# INLINE filterCrossHighLevel #-} filterCrossHighLevel :: (PrimMonad m, Functor m, Num a, Mult a b, VG.Vector v a, VG.Vector v b, VGM.MVector vm a) => v b -> Int -> v a -> v a -> vm (PrimState m) a -> m ()-filterCrossHighLevel coeffs num lastBuf nextBuf outBuf = fill (VFSM.generate num dotProd) outBuf+filterCrossHighLevel coeffs num lastBuf nextBuf outBuf = fill (VFB.generate num dotProd) outBuf     where     dotProd i = VG.sum $ VG.zipWith mult (VG.unsafeDrop i lastBuf VG.++ nextBuf) coeffs 
hs_sources/SDR/PipeUtils.hs view
@@ -54,7 +54,11 @@             rate' (buffers + 1)     rate' 1 -pMapAccum :: (Monad m) => (acc -> x -> (acc, y)) -> acc -> Pipe x y m () +-- | mapAccum for Pipes+pMapAccum :: (Monad m) +          => (acc -> x -> (acc, y)) -- ^ Accumulating function+          -> acc                    -- ^ Initial value of the accumulator+          -> Pipe x y m ()  pMapAccum func acc = go acc     where     go acc = do
hs_sources/SDR/Serialize.hs view
@@ -31,8 +31,8 @@ import           Data.ByteString          as BS import           System.IO -import           Data.Vector.Generic      as VG hiding ((++))-import           Data.Vector.Storable     as VS hiding ((++))+import qualified Data.Vector.Generic      as VG hiding ((++))+import qualified Data.Vector.Storable     as VS hiding ((++))  import           Pipes import qualified Pipes.Prelude            as P
hs_sources/SDR/Util.hs view
@@ -6,42 +6,72 @@     Mult,     mult, -    -- * Conversion to Floating Point+    -- * Conversion to floating point for reception+    -- ** RTLSDR     interleavedIQUnsigned256ToFloat,     interleavedIQUnsignedByteToFloat,     interleavedIQUnsignedByteToFloatSSE,     interleavedIQUnsignedByteToFloatAVX,     interleavedIQUnsignedByteToFloatFast, +    -- ** BladeRF     interleavedIQSigned2048ToFloat,     interleavedIQSignedWordToFloat,     interleavedIQSignedWordToFloatSSE,     interleavedIQSignedWordToFloatAVX,+    interleavedIQSignedWordToFloatFast, +    -- * Conversion from floating point for transmission+    -- ** BladeRF+    complexFloatToInterleavedIQSigned2048,+    complexFloatToInterleavedIQSignedWord,+     -- * Scaling     scaleC,     scaleCSSE,     scaleCAVX,     scaleFast, -    -- * Misc Utils+    -- * Mapping over complex numbers     cplxMap,-    quarterBandUp++    -- * Frequency shifting+    halfBandUp,+    quarterBandUp,++    -- * Data streams+    streamString,+    streamRandom,++    -- * Automatic gain control+    agc,+    agcPipe     ) where  import           Foreign.C.Types import           Data.Complex-import           Data.Vector.Generic          as VG   hiding ((++))+import qualified Data.Vector.Generic          as VG    import qualified Data.Vector.Generic.Mutable  as VGM-import           Data.Vector.Storable         as VS   hiding ((++))-import           Data.Vector.Storable.Mutable as VSM  +import qualified Data.Vector.Storable         as VS   +import qualified Data.Vector.Storable.Mutable as VSM   import           Control.Monad.Primitive import           Unsafe.Coerce import           Foreign.Ptr import           System.IO.Unsafe import           Foreign.Storable.Complex+import           Control.Monad+import qualified System.Random.MWC as R+import           Data.Bits+import           Pipes+import qualified Pipes.Prelude as P+import           Data.Word+import           Foreign.Storable+import           Control.Arrow as A+import           Data.Tuple  import           SDR.CPUID+import           SDR.VectorUtils+import           SDR.PipeUtils  -- | A class for things that can be multiplied by a scalar. class Mult a b where@@ -116,7 +146,7 @@ foreign import ccall unsafe "convertCBladeRF"     convertCBladeRF_c :: CInt -> Ptr CShort -> Ptr CFloat -> IO () --- | Same as `interleavedIQUnsigned256ToFloat` but written in C and specialized for unsigned byte inputs and Float outputs.+-- | Same as `interleavedIQUnsigned256ToFloat` but written in C and specialized for signed short inputs and Float outputs. interleavedIQSignedWordToFloat :: VS.Vector CShort -> VS.Vector (Complex Float) interleavedIQSignedWordToFloat inBuf = unsafePerformIO $ do     outBuf <- VGM.new $ VG.length inBuf `quot` 2@@ -128,7 +158,7 @@ foreign import ccall unsafe "convertCSSEBladeRF"     convertCSSEBladeRF_c :: CInt -> Ptr CShort -> Ptr CFloat -> IO () --- | Same as `interleavedIQUnsigned256ToFloat` but written in C using SSE intrinsics and specialized for unsigned byte inputs and Float outputs.+-- | Same as `interleavedIQUnsigned256ToFloat` but written in C using SSE intrinsics and specialized for signed short inputs and Float outputs. interleavedIQSignedWordToFloatSSE :: VS.Vector CShort -> VS.Vector (Complex Float) interleavedIQSignedWordToFloatSSE inBuf = unsafePerformIO $ do     outBuf <- VGM.new $ VG.length inBuf `quot` 2@@ -140,7 +170,7 @@ foreign import ccall unsafe "convertCAVXBladeRF"     convertCAVXBladeRF_c :: CInt -> Ptr CShort -> Ptr CFloat -> IO () --- | Same as `interleavedIQUnsigned256ToFloat` but written in C using AVX intrinsics and specialized for unsigned byte inputs and Float outputs.+-- | Same as `interleavedIQUnsigned256ToFloat` but written in C using AVX intrinsics and specialized for signed short inputs and Float outputs. interleavedIQSignedWordToFloatAVX :: VS.Vector CShort -> VS.Vector (Complex Float) interleavedIQSignedWordToFloatAVX inBuf = unsafePerformIO $ do     outBuf <- VGM.new $ VG.length inBuf `quot` 2@@ -149,6 +179,33 @@             convertCAVXBladeRF_c (fromIntegral $ VG.length inBuf) iPtr oPtr     VG.freeze outBuf +-- | Same as `interleavedIQSigned2048ToFloat` but uses the fastest SIMD instruction set your processor supports and specialized for signed short inputs and Float outputs.+interleavedIQSignedWordToFloatFast :: CPUInfo -> VS.Vector CShort -> VS.Vector (Complex Float)+interleavedIQSignedWordToFloatFast info = featureSelect info interleavedIQSignedWordToFloat [(hasAVX2, interleavedIQSignedWordToFloatAVX), (hasSSE42, interleavedIQSignedWordToFloatSSE)]++-- | Create a vector of interleaved I Q component integral samples from a vector of complex Floats. Each input ranges from -2048 to 2047. This is the format the BladeRF uses.+complexFloatToInterleavedIQSigned2048 :: (Integral b, RealFrac a, VG.Vector v1 (Complex a), VG.Vector v2 b) => v1 (Complex a) -> v2 b+complexFloatToInterleavedIQSigned2048 input = VG.generate (VG.length input * 2) convert+    where+    {-# INLINE convert #-}+    convert idx  +        | even idx = convert' $ realPart (input `VG.unsafeIndex` (idx `quot` 2))+        | odd  idx = convert' $ imagPart (input `VG.unsafeIndex` (idx `quot` 2))+    {-# INLINE convert' #-}+    convert' val = round $ val * 2048++foreign import ccall unsafe "convertBladeRFTransmit"+    convertBladeRFTransmit_c :: CInt -> Ptr CFloat -> Ptr CShort -> IO ()++-- | Same as `complexFloatToInterleavedIQSigned2048` but written in C and specialized for Float inputs and signed short outputs.+complexFloatToInterleavedIQSignedWord :: VS.Vector (Complex Float) -> VS.Vector CShort+complexFloatToInterleavedIQSignedWord inBuf = unsafePerformIO $ do+    outBuf <- VGM.new $ VG.length inBuf * 2+    VS.unsafeWith (unsafeCoerce inBuf) $ \iPtr -> +        VSM.unsafeWith outBuf $ \oPtr -> +            convertBladeRFTransmit_c (fromIntegral $ VG.length inBuf * 2) iPtr oPtr+    VG.freeze outBuf+     -- | Scaling foreign import ccall unsafe "scale"     scale_c :: CInt -> CFloat -> Ptr CFloat -> Ptr CFloat -> IO ()@@ -199,6 +256,16 @@         -> Complex b -- ^ Output complex number cplxMap f (x :+ y) = f x :+ f y +-- | Multiplication by this vector shifts all frequencies up by 1/2 of the sampling frequency+halfBandUp :: (VG.Vector v n, Num n) +           => Int -- ^ The length of the Vector+           -> v n +halfBandUp size = VG.generate size func+    where+    func idx +        | even idx  = 1+        | otherwise = -1+ -- | Multiplication by this vector shifts all frequencies up by 1/4 of the sampling frequency quarterBandUp :: (VG.Vector v (Complex n), Num n)                => Int -- ^ The length of the Vector@@ -212,4 +279,67 @@         | m == 3 = 0    :+ (-1)         where         m = idx `mod` 4++-- | A Producer that streams vectors of the bits that make up the string argument concatenated repeatedly. Each bit is encoded as a float with value (+1) for 1 and (-1) for 0.+streamString :: forall m b. (FiniteBits b, Monad m) +             => [b] -- ^ The string whose bits are to be streamed+             -> Int -- ^ The size of each streamed vector+             -> Producer (VS.Vector Float) m ()+streamString str size = P.unfoldr (return . Right . func) (str, 0)+    where+    bitsPerChar = finiteBitSize (undefined :: b)+    toFloat :: Bool -> Float+    toFloat x   = if x then 1 else (-1)+    func = vUnfoldr size funcy +        where+        funcy ([], offsetChar) = funcy (str, 0)+        funcy (rem@(x:xs), offsetChar)+            | offsetChar == bitsPerChar = funcy (xs, 0)+            | otherwise                 = (toFloat $ testBit x offsetChar, (rem, offsetChar + 1))++-- | A Producer that streams vectors of random bits. Each bit is encoded as a float with value (+1) for 1 and (-1) for 0.+streamRandom :: forall m. PrimMonad m +             => Int -- ^ The size of each streamed vector+             -> Producer (VS.Vector Float) m ()+streamRandom size = do+    gen   <- lift R.create +    start <- lift $ R.uniform gen+    P.unfoldr (liftM Right . func gen) (start, 0)+    where+    toFloat :: Bool -> Float+    toFloat x   = if x then 1 else (-1)+    func :: R.Gen (PrimState m) -> (Word64, Int) -> m (VS.Vector Float, (Word64, Int))+    func gen = vUnfoldrM size funcy +        where+        funcy (current, offset) = do+            let res =  toFloat $ testBit current offset+            if offset == 63 then do+                current' <- R.uniform gen+                return (res, (current', 0))+            else return (res, (current, offset+1))++(a :+ b) `cdiv` y = (a/y) :+ (b/y)+(a :+ b) `cmul` y = (a*y) :+ (b*y)++-- | Simple automatic gain control +agc :: (Num a, Storable a, RealFloat a) +    => a                          -- ^ a+    -> a                          -- ^ reference+    -> a                          -- ^ initial state+    -> VS.Vector (Complex a)      -- ^ input vector+    -> (a, VS.Vector (Complex a)) -- ^ (final state, output vector)+agc mu reference state input = A.first snd $ swap $ vUnfoldr (VS.length input) go (0, state)+    where+    go (offset, state) = +        let+            corrected = (input VS.! offset) `cmul` state+            state'    = state + mu * (reference - magnitude corrected)+        in  (corrected, (offset + 1, state'))++-- | Simple automatic gain control pipe+agcPipe :: (Num a, Storable a, RealFloat a, Monad m)+        => a -- ^ a+        -> a -- ^ reference+        -> Pipe (VS.Vector (Complex a)) (VS.Vector (Complex a)) m ()+agcPipe mu reference = pMapAccum (agc mu reference) 1 
hs_sources/SDR/VectorUtils.hs view
@@ -5,18 +5,19 @@     mapAccumMV,     stride,     fill,-    vUnfoldr+    copyInto,+    vUnfoldr,+    vUnfoldrM     ) where  import Control.Monad import Control.Monad.Primitive import Control.Monad.ST -import           Data.Vector.Generic               as VG   hiding ((++))+import qualified Data.Vector.Generic               as VG    import qualified Data.Vector.Generic.Mutable       as VGM-import           Data.Vector.Fusion.Stream.Monadic         hiding ((++))-import qualified Data.Vector.Fusion.Stream         as VFS  hiding ((++))-import qualified Data.Vector.Fusion.Stream.Monadic as VFSM hiding ((++))+import qualified Data.Vector.Fusion.Bundle         as VFB  +import           Data.Vector.Fusion.Stream.Monadic as VFSM  {-| Like mapAccumL but monadic and over vectors. Doesn't return the     accumulator at the end because it doesn't seem to be possible to do@@ -25,18 +26,18 @@ mapAccumMV :: (Monad m)             => (acc -> x -> m (acc, y)) -- ^ The function            -> acc                      -- ^ The initial accumulator-           -> Stream m x               -- ^ The input stream+           -> VFSM.Stream m x          -- ^ The input stream            -> Stream m y               -- ^ The output stream-mapAccumMV func z (Stream step s sz) = Stream step' (s, z) sz+mapAccumMV func z (VFSM.Stream step s) = VFSM.Stream step' (s, z)     where     step' (s, acc) = do         r <- step s         case r of-            Yield y s' -> do+            VFB.Yield y s' -> do                 (!acc', !res) <- func acc y -                return $ Yield res (s', acc')-            Skip    s' -> return $ Skip (s', acc)-            Done       -> return Done+                return $ VFB.Yield res (s', acc')+            VFB.Skip    s' -> return $ VFB.Skip (s', acc)+            VFB.Done       -> return VFB.Done  {-| Create a vector from another vector containing only the elements that     occur every stride elements in the source vector.@@ -46,25 +47,39 @@        => Int -- ^ The stride        -> v a -- ^ The input Vector        -> v a -- ^ The output Vector-stride str inv = VG.unstream $ VFS.unfoldr func 0+stride str inv = VG.unstream $ VFB.unfoldr func 0     where     len = VG.length inv     func i | i >= len  = Nothing            | otherwise = Just (VG.unsafeIndex inv i, i + str) --- | Fill a mutable vector from a monadic stream+-- | Fill a mutable vector from a monadic stream. This appears to be missing from the Vector library. {-# INLINE fill #-} fill :: (PrimMonad m, Functor m, VGM.MVector vm a) -     => VFS.MStream m a    -- ^ The input Stream+     => VFB.Bundle v a    -- ^ The input Stream      -> vm (PrimState m) a -- ^ The mutable Vector to stream into      -> m ()-fill str outBuf = void $ VFSM.foldM' put 0 str+fill str outBuf = void $ VFB.foldM' put 0 str     where      put i x = do         VGM.unsafeWrite outBuf i x         return $ i + 1 -vUnfoldr :: VG.Vector v x => Int -> (acc -> (x, acc)) -> acc -> (v x, acc)+-- | Copy a Vector into a mutable vector+{-# INLINE copyInto #-}+copyInto :: (PrimMonad m, VGM.MVector vm a, VG.Vector v a)+         => vm (PrimState m) a -- ^ The destination+         -> v a                -- ^ The source+         -> m ()+copyInto dst src = fill (VG.stream src) dst++-- | Similar to unfoldrN from the vector package but the generator function cannot terminate and it returns the final value of the seed in addition to the vector+{-# INLINE vUnfoldr #-}+vUnfoldr :: VG.Vector v x +         => Int               -- ^ Generates a vector with this size+         -> (acc -> (x, acc)) -- ^ The generator function+         -> acc               -- ^ The initial value of the seed+         -> (v x, acc)        -- ^ The (vector, final value of seed) result vUnfoldr size func acc = runST $ do     vect <- VGM.new size     acc' <- go vect 0 acc@@ -79,3 +94,26 @@                 let (res, acc') = func acc                 VGM.write vect offset res                 go' (offset + 1) acc'++-- | The same as `vUnfoldr` but the generator function is monadic+{-# INLINE vUnfoldrM #-}+vUnfoldrM :: (PrimMonad m, VG.Vector v x) +          => Int                 -- ^ Generates a vector with this size+          -> (acc -> m (x, acc)) -- ^ The monadic generator function                    +          -> acc                 -- ^ The initial value of the seed+          -> m (v x, acc)        -- ^ The (vector, final value of seed) result+vUnfoldrM size func acc = do+    vect <- VGM.new size+    acc' <- go vect 0 acc+    vect' <- VG.unsafeFreeze vect+    return (vect', acc')+    where+    go vect offset acc = go' offset acc+        where+        go' offset acc +            | offset == size = return acc+            | otherwise      = do+                (res, acc') <- func acc+                VGM.write vect offset res+                go' (offset + 1) acc'+
sdr.cabal view
@@ -1,5 +1,5 @@ name:                sdr-version:             0.1.0.5+version:             0.1.0.6 synopsis:            A software defined radio library description:              Write software defined radio applications in Haskell.@@ -24,7 +24,7 @@     .     * PulseAudio sound sink     .-    * <http://sdr.osmocom.org/trac/wiki/rtl-sdr rtl-sdr> based radio source supported and other sources are easily added+    * <http://sdr.osmocom.org/trac/wiki/rtl-sdr rtl-sdr> and <https://nuand.com BladeRF> based radio sources/sinks supported and other sources are easily added     .     See <https://github.com/adamwalker/sdr> for more features and screenshots.     .@@ -80,7 +80,7 @@         pipes-bytestring     >=2.0   && <2.2,          dynamic-graph        ==0.1.0.8,         array                >=0.4   && <0.6, -        vector               >=0.10  && <0.11,+        vector               >=0.11  && <0.12,         tuple                >=0.2   && <0.4,          OpenGL               >=2.11  && <2.13,          GLFW-b               >=1.4.7 && <1.4.8,@@ -89,11 +89,12 @@         pango                >=0.13  && <0.14,          containers           >=0.5   && <0.6,          cairo                >=0.13  && <0.14, -        cereal               >=0.4   && <0.5, -        optparse-applicative >=0.11  && <0.12, +        cereal               >=0.4   && <0.6, +        optparse-applicative >=0.11  && <0.13,          Decimal              >=0.4   && <0.5, -        Chart                >=1.3   && <1.5, -        Chart-cairo          >=1.3   && <1.5+        Chart                >=1.3   && <1.6, +        Chart-cairo          >=1.3   && <1.6,+        mwc-random     -- hs-source-dirs:           default-language:    Haskell2010     ghc-options:         -O2@@ -115,8 +116,8 @@     build-depends:                base                       >=4.6  && <4.9,          QuickCheck                 >=2.8  && <2.9, -        vector                     >=0.10 && <0.11, -        sdr                        ==0.1.0.5, +        vector                     >=0.11 && <0.12, +        sdr                        ==0.1.0.6,          primitive                  >=0.5  && <0.7,          storable-complex           >=0.2  && <0.3,         test-framework             >=0.8  && <0.9,@@ -131,8 +132,8 @@     build-depends:                base             >=4.6  && <4.9,          criterion        >=1.0  && <1.2,-        vector           >=0.10 && <0.11, -        sdr              ==0.1.0.5, +        vector           >=0.11 && <0.12, +        sdr              ==0.1.0.6,          primitive        >=0.5  && <0.7,          storable-complex >=0.2  && <0.3     hs-source-dirs:      benchmarks
tests/TestSuite.hs view
@@ -39,7 +39,8 @@             testProperty "complex" propDecimationComplex         ],         testGroup "resamplers" [-            testProperty "real" propResamplingReal+            testProperty "real" propResamplingReal,+            testProperty "complex" propResamplingComplex         ],         testGroup "conversion" [             testProperty "rtlsdr"  propConversionRTLSDR,@@ -188,6 +189,37 @@         res <- run $ sameResultM eqDelta $ map (getResult num $) $ hasFeatures [                 (const True, void . resampleHighLevel       interpolation decimation vCoeffs offset num vInput),                 (const True, void . resampleCRR             num interpolation decimation offset vCoeffs vInput),+                (const True, void . resampler3              group num vInput),+                (hasSSE42,   void . resampler4              group num vInput),+                (hasAVX,     void . resampler5              group num vInput)+            ]+        assert res++    propResamplingComplex = +        forAll sizes $ \size -> +            forAll (vectorOf size (choose (-10, 10))) $ \inBufR -> +                forAll (vectorOf size (choose (-10, 10))) $ \inBufI -> +                    forAll (elements [32 .. 512]) $ \numCoeffs -> +                        forAll (vectorOf numCoeffs (choose (-10, 10))) $ \coeffs -> +                           forAll (elements $ tail factors') $ \decimation -> +                               forAll (elements $ filter (< decimation) factors') $ \interpolation -> +                                   forAll (elements [0..interpolation - 1]) $ \group -> +                                       testResamplingComplex size group numCoeffs interpolation decimation coeffs $ zipWith (:+) inBufR inBufI++    testResamplingComplex :: Int -> Int -> Int -> Int -> Int -> [Float] -> [Complex Float] -> Property+    testResamplingComplex size group numCoeffs interpolation decimation coeffs inBuf = monadicIO $ do+        let vCoeffsHalf = VS.fromList coeffs+            vCoeffs     = VS.fromList $ coeffs ++ reverse coeffs+            vInput      = VS.fromList inBuf+            num         = (size - numCoeffs*2 + 1) `quot` decimation+            offset       = interpolation - 1 - ((interpolation + group * decimation - 1) `mod` interpolation) ++        resampler3 <- run $ resampleCRC    interpolation decimation (coeffs ++ reverse coeffs)+        resampler4 <- run $ resampleCSSERC interpolation decimation (coeffs ++ reverse coeffs)+        resampler5 <- run $ resampleCAVXRC interpolation decimation (coeffs ++ reverse coeffs)++        res <- run $ sameResultM eqDeltaC $ map (getResult num $) $ hasFeatures [+                (const True, void . resampleHighLevel       interpolation decimation vCoeffs offset num vInput),                 (const True, void . resampler3              group num vInput),                 (hasSSE42,   void . resampler4              group num vInput),                 (hasAVX,     void . resampler5              group num vInput)