fec (empty) → 0.1
raw patch · 6 files changed
+1197/−0 lines, 6 filesdep +basedep +bytestringsetup-changed
Dependencies added: base, bytestring
Files
- README.txt +239/−0
- Setup.lhs +3/−0
- fec.cabal +30/−0
- haskell/Codec/FEC.hs +246/−0
- zfec/fec.c +570/−0
- zfec/fec.h +109/−0
+ README.txt view
@@ -0,0 +1,239 @@+ * Intro and Licence++This package implements an "erasure code", or "forward error correction+code".++You may use this package under the GNU General Public License, version 2 or, at+your option, any later version. You may use this package under the Transitive+Grace Period Public Licence, version 1.0. (You may choose to use this package+under the terms of either licence, at your option.) See the file COPYING.GPL+for the terms of the GNU General Public License, version 2. See the file+COPYING.TGPPL.html for the terms of the Transitive Grace Period Public Licence,+version 1.0.++The most widely known example of an erasure code is the RAID-5 algorithm+which makes it so that in the event of the loss of any one hard drive, the+stored data can be completely recovered. The algorithm in the zfec package+has a similar effect, but instead of recovering from the loss of only a+single element, it can be parameterized to choose in advance the number of+elements whose loss it can tolerate.++This package is largely based on the old "fec" library by Luigi Rizzo et al.,+which is a mature and optimized implementation of erasure coding. The zfec+package makes several changes from the original "fec" package, including+addition of the Python API, refactoring of the C API to support zero-copy+operation, a few clean-ups and optimizations of the core code itself, and the+addition of a command-line tool named "zfec".+++ * Installation++This package is managed with the "setuptools" package management tool. To+build and install the package directly into your system, just run "python+./setup.py install". If you prefer to keep the package limited to a specific+directory so that you can manage it yourself (perhaps by using the "GNU+stow") tool, then give it these arguments: "python ./setup.py install+--single-version-externally-managed+--record=${specificdirectory}/zfec-install.log --prefix=${specificdirectory}"++To run the self-tests, execute "python ./setup.py test" (or if you have +Twisted Python installed, you can run "trial zfec" for nicer output and test +options.)+++ * Community++The source is currently available via darcs on the web with the command:++darcs get http://allmydata.org/source/zfec++More information on darcs is available at http://darcs.net++Please join the zfec mailing list and submit patches:++<http://allmydata.org/cgi-bin/mailman/listinfo/zfec-dev>+++ * Overview++This package performs two operations, encoding and decoding. Encoding takes+some input data and expands its size by producing extra "check blocks", also+called "secondary blocks". Decoding takes some data -- any combination of+blocks of the original data (called "primary blocks") and "secondary blocks",+and produces the original data.++The encoding is parameterized by two integers, k and m. m is the total number+of blocks produced, and k is how many of those blocks are necessary to+reconstruct the original data. m is required to be at least 1 and at most 256,+and k is required to be at least 1 and at most m.++(Note that when k == m then there is no point in doing erasure coding -- it+degenerates to the equivalent of the Unix "split" utility which simply splits+the input into successive segments. Similarly, when k == 1 it degenerates to+the equivalent of the unix "cp" utility -- each block is a complete copy of the+input data. The "zfec" command-line tool does not implement these degenerate +cases.)++Note that each "primary block" is a segment of the original data, so its size+is 1/k'th of the size of original data, and each "secondary block" is of the+same size, so the total space used by all the blocks is m/k times the size of+the original data (plus some padding to fill out the last primary block to be+the same size as all the others). In addition to the data contained in the +blocks themselves there are also a few pieces of metadata which are necessary +for later reconstruction. Those pieces are: 1. the value of K, 2. the value +of M, 3. the sharenum of each block, 4. the number of bytes of padding +that were used. The "zfec" command-line tool compresses these pieces of data +and prepends them to the beginning of each share, so each the sharefile +produced by the "zfec" command-line tool is between one and four bytes larger +than the share data alone.++The decoding step requires as input k of the blocks which were produced by the+encoding step. The decoding step produces as output the data that was earlier+input to the encoding step.+++ * Command-Line Tool++NOTE: the format of the sharefiles was changed in zfec v1.1 to allow K == 1 +and K == M. This change of the format of sharefiles means that zfec >= v1.1 +cannot read sharefiles produced by zfec < v1.1.++The bin/ directory contains two Unix-style, command-line tools "zfec" and +"zunfec". Execute "zfec --help" or "zunfec --help" for usage instructions.++Note: a Unix-style tool like "zfec" does only one thing -- in this case+erasure coding -- and leaves other tasks to other tools. Other Unix-style+tools that go well with zfec include "GNU tar" for archiving multiple files+and directories into one file, "rzip" or "lrzip" for compression, and "GNU+Privacy Guard" for encryption or "sha256sum" for integrity. It is important+to do things in order: first archive, then compress, then either encrypt or+sha256sum, then erasure code. Note that if GNU Privacy Guard is used for+privacy, then it will also ensure integrity, so the use of sha256sum is+unnecessary in that case.+++ * Performance Measurements++On my Athlon 64 2.4 GHz workstation (running Linux), the "zfec" command-line+tool encoded a 160 MB file with m=100, k=94 (about 6% redundancy) in 3.9+seconds, where the "par2" tool encoded the file with about 6% redundancy in+27 seconds. zfec encoded the same file with m=12, k=6 (100% redundancy) in+4.1 seconds, where par2 encoded it with about 100% redundancy in 7 minutes+and 56 seconds.++The underlying C library in benchmark mode encoded from a file at about +4.9 million bytes per second and decoded at about 5.8 million bytes per second.++On Peter's fancy Intel Mac laptop (2.16 GHz Core Duo), it encoded from a file+at about 6.2 million bytes per second.++On my even fancier Intel Mac laptop (2.33 GHz Core Duo), it encoded from a file+at about 6.8 million bytes per second.++On my old PowerPC G4 867 MHz Mac laptop, it encoded from a file at about 1.3+million bytes per second.+++ * API++Each block is associated with "blocknum". The blocknum of each primary block is+its index (starting from zero), so the 0'th block is the first primary block,+which is the first few bytes of the file, the 1'st block is the next primary+block, which is the next few bytes of the file, and so on. The last primary+block has blocknum k-1. The blocknum of each secondary block is an arbitrary+integer between k and 255 inclusive. (When using the Python API, if you don't+specify which blocknums you want for your secondary blocks when invoking+encode(), then it will by default provide the blocks with ids from k to m-1+inclusive.)++ ** C API++fec_encode() takes as input an array of k pointers, where each pointer points+to a memory buffer containing the input data (i.e., the i'th buffer contains+the i'th primary block). There is also a second parameter which is an array of+the blocknums of the secondary blocks which are to be produced. (Each element+in that array is required to be the blocknum of a secondary block, i.e. it is+required to be >= k and < m.)++The output from fec_encode() is the requested set of secondary blocks which are+written into output buffers provided by the caller.++fec_decode() takes as input an array of k pointers, where each pointer points+to a buffer containing a block. There is also a separate input parameter which+is an array of blocknums, indicating the blocknum of each of the blocks which is+being passed in.++The output from fec_decode() is the set of primary blocks which were missing+from the input and had to be reconstructed. These reconstructed blocks are+written into output buffers provided by the caller.++ ** Python API++encode() and decode() take as input a sequence of k buffers, where a "sequence"+is any object that implements the Python sequence protocol (such as a list or+tuple) and a "buffer" is any object that implements the Python buffer protocol+(such as a string or array). The contents that are required to be present in+these buffers are the same as for the C API.++encode() also takes a list of desired blocknums. Unlike the C API, the Python+API accepts blocknums of primary blocks as well as secondary blocks in its list+of desired blocknums. encode() returns a list of buffer objects which contain+the blocks requested. For each requested block which is a primary block, the+resulting list contains a reference to the apppropriate primary block from the+input list. For each requested block which is a secondary block, the list+contains a newly created string object containing that block.++decode() also takes a list of integers indicating the blocknums of the blocks+being passed int. decode() returns a list of buffer objects which contain all+of the primary blocks of the original data (in order). For each primary block+which was present in the input list, then the result list simply contains a+reference to the object that was passed in the input list. For each primary+block which was not present in the input, the result list contains a newly+created string object containing that primary block.++Beware of a "gotcha" that can result from the combination of mutable data and+the fact that the Python API returns references to inputs when possible.++Returning references to its inputs is efficient since it avoids making an+unnecessary copy of the data, but if the object which was passed as input is+mutable and if that object is mutated after the call to zfec returns, then the+result from zfec -- which is just a reference to that same object -- will also+be mutated. This subtlety is the price you pay for avoiding data copying. If+you don't want to have to worry about this then you can simply use immutable+objects (e.g. Python strings) to hold the data that you pass to zfec.+++ * Utilities++The filefec.py module has a utility function for efficiently reading a file+and encoding it piece by piece. This module is used by the "zfec" and +"zunfec" command-line tools from the bin/ directory.+++ * Dependencies++A C compiler is required. To use the Python API or the command-line tools a+Python interpreter is also required. We have tested it with Python v2.4 and+v2.5.+++ * Acknowledgements++Thanks to the author of the original fec lib, Luigi Rizzo, and the folks that+contributed to it: Phil Karn, Robert Morelos-Zaragoza, Hari Thirumoorthy, and+Dan Rubenstein. Thanks to the Mnet hackers who wrote an earlier Python+wrapper, especially Myers Carpenter and Hauke Johannknecht. Thanks to Brian+Warner and Amber O'Whielacronx for help with the API, documentation, +debugging, compression, and unit tests. Thanks to the creators of GCC +(starting with Richard M. Stallman) and Valgrind (starting with Julian Seward) +for a pair of excellent tools. Thanks to my coworkers at Allmydata -- +http://allmydata.com -- Fabrice Grinda, Peter Secor, Rob Kinninmont, Brian +Warner, Zandr Milewski, Justin Boreta, Mark Meras for sponsoring this work and +releasing it under a Free Software licence.+++Enjoy!++Zooko Wilcox-O'Hearn+2007-10-01+Boulder, Colorado
+ Setup.lhs view
@@ -0,0 +1,3 @@+#!/usr/bin/env runhaskell+> import Distribution.Simple+> main = defaultMain
+ fec.cabal view
@@ -0,0 +1,30 @@+name: fec+version: 0.1+license: GPL+license-file: README.txt+author: Adam Langley <agl@imperialviolet.org>+maintainer: Adam Langley <agl@imperialviolet.org>+description: This code, based on zfec by Zooko, based on code by Luigi+ Rizzo implements an erasure code, or forward error+ correction code. The most widely known example of an erasure+ code is the RAID-5 algorithm which makes it so that in the+ event of the loss of any one hard drive, the stored data can+ be completely recovered. The algorithm in the zfec package+ has a similar effect, but instead of recovering from the loss+ of only a single element, it can be parameterized to choose in+ advance the number of elements whose loss it can tolerate.+build-type: Simple+homepage: http://allmydata.org/source/zfec+synopsis: Forward error correction of ByteStrings+category: Codec+build-depends: base, bytestring>=0.9+stability: provisional+tested-with: GHC == 6.8.2+exposed-modules: Codec.FEC+extensions: ForeignFunctionInterface+hs-source-dirs: haskell+ghc-options: -Wall+c-sources: zfec/fec.c+cc-options: -std=c99+include-dirs: zfec+extra-source-files: zfec/fec.h
+ haskell/Codec/FEC.hs view
@@ -0,0 +1,246 @@+{-# LANGUAGE ForeignFunctionInterface, EmptyDataDecls #-}+-- |+-- Module: Codec.FEC+-- Copyright: Adam Langley+-- License: BSD3+--+-- Stability: experimental+--+-- The module provides k of n encoding - a way to generate (n - k) secondary+-- blocks of data from k primary blocks such that any k blocks (primary or+-- secondary) are sufficient to regenerate all blocks.+--+-- All blocks must be the same length and you need to keep track of which+-- blocks you have in order to tell decode. By convention, the blocks are+-- numbered 0..(n - 1) and blocks numbered < k are the primary blocks.++module Codec.FEC (+ FECParams+ , fec+ , encode+ , decode++ -- * Utility functions+ , secureDivide+ , secureCombine+ , enFEC+ , deFEC+ ) where++import qualified Data.ByteString as B+import qualified Data.ByteString.Unsafe as BU+import qualified Data.ByteString.Internal as BI+import Data.Word (Word8)+import Data.Bits (xor)+import Data.List (sortBy, partition, (\\), nub)+import Foreign.Ptr+import Foreign.Storable (sizeOf, poke)+import Foreign.ForeignPtr+import Foreign.C.Types+import Foreign.Marshal.Alloc+import Foreign.Marshal.Array (withArray, advancePtr)+import System.IO (withFile, IOMode(..))+import System.IO.Unsafe (unsafePerformIO)++data CFEC+data FECParams = FECParams (ForeignPtr CFEC) Int Int++instance Show FECParams where+ show (FECParams _ k n) = "FEC (" ++ show k ++ ", " ++ show n ++ ")"++foreign import ccall unsafe "fec_new" _new :: CUInt -- ^ k+ -> CUInt -- ^ n+ -> IO (Ptr CFEC)+foreign import ccall unsafe "&fec_free" _free :: FunPtr (Ptr CFEC -> IO ())+foreign import ccall unsafe "fec_encode" _encode :: Ptr CFEC+ -> Ptr (Ptr Word8) -- ^ primary blocks+ -> Ptr (Ptr Word8) -- ^ (output) secondary blocks+ -> Ptr CUInt -- ^ array of secondary block ids+ -> CSize -- ^ length of previous+ -> CSize -- ^ block length+ -> IO ()+foreign import ccall unsafe "fec_decode" _decode :: Ptr CFEC+ -> Ptr (Ptr Word8) -- ^ input blocks+ -> Ptr (Ptr Word8) -- ^ output blocks+ -> Ptr CUInt -- ^ array of input indexes+ -> CSize -- ^ block length+ -> IO ()++-- | Return true if the given @k@ and @n@ values are valid+isValidConfig :: Int -> Int -> Bool+isValidConfig k n+ | k >= n = False+ | k < 1 = False+ | n < 1 = False+ | n > 255 = False+ | otherwise = True++-- | Return a FEC with the given parameters.+fec :: Int -- ^ the number of primary blocks+ -> Int -- ^ the total number blocks, must be < 256+ -> FECParams+fec k n =+ if not (isValidConfig k n)+ then error $ "Invalid FEC parameters: " ++ show k ++ " " ++ show n+ else unsafePerformIO (do+ cfec <- _new (fromIntegral k) (fromIntegral n)+ params <- newForeignPtr _free cfec+ return $ FECParams params k n)++-- | Create a C array of unsigned from an input array+uintCArray :: [Int] -> ((Ptr CUInt) -> IO a) -> IO a+uintCArray xs f = withArray (map fromIntegral xs) f++-- | Convert a list of ByteStrings to an array of pointers to their data+byteStringsToArray :: [B.ByteString] -> ((Ptr (Ptr Word8)) -> IO a) -> IO a+byteStringsToArray inputs f = do+ let l = length inputs+ allocaBytes (l * sizeOf (undefined :: Ptr Word8)) (\array -> do+ let inner _ [] = f array+ inner array' (bs : bss) = BU.unsafeUseAsCString bs (\ptr -> do+ poke array' $ castPtr ptr+ inner (advancePtr array' 1) bss)+ inner array inputs)++-- | Return True iff all the given ByteStrings are the same length+allByteStringsSameLength :: [B.ByteString] -> Bool+allByteStringsSameLength [] = True+allByteStringsSameLength (bs : bss) = all ((==) (B.length bs)) $ map B.length bss++-- | Run the given function with a pointer to an array of @n@ pointers to+-- buffers of size @size@. Return these buffers as a list of ByteStrings+createByteStringArray :: Int -- ^ the number of buffers requested+ -> Int -- ^ the size of each buffer+ -> ((Ptr (Ptr Word8)) -> IO ())+ -> IO [B.ByteString]+createByteStringArray n size f = do+ allocaBytes (n * sizeOf (undefined :: Ptr Word8)) (\array -> do+ allocaBytes (n * size) (\ptr -> do+ mapM_ (\i -> poke (advancePtr array i) (advancePtr ptr (size * i))) [0..(n - 1)]+ f array+ mapM (\i -> B.packCStringLen (castPtr $ advancePtr ptr (i * size), size)) [0..(n - 1)]))++-- | Generate the secondary blocks from a list of the primary blocks. The+-- primary blocks must be in order and all of the same size. There must be+-- @k@ primary blocks.+encode :: FECParams+ -> [B.ByteString] -- ^ a list of @k@ input blocks+ -> [B.ByteString] -- ^ (n - k) output blocks+encode (FECParams params k n) inblocks+ | length inblocks /= k = error "Wrong number of blocks to FEC encode"+ | not (allByteStringsSameLength inblocks) = error "Not all inputs to FEC encode are the same length"+ | otherwise = unsafePerformIO (do+ let sz = B.length $ head inblocks+ withForeignPtr params (\cfec -> do+ byteStringsToArray inblocks (\src -> do+ createByteStringArray (n - k) sz (\fecs -> do+ uintCArray [k..(n - 1)] (\block_nums -> do+ _encode cfec src fecs block_nums (fromIntegral (n - k)) $ fromIntegral sz)))))++-- | A sort function for tagged assoc lists+sortTagged :: [(Int, a)] -> [(Int, a)]+sortTagged = sortBy (\a b -> compare (fst a) (fst b))++-- | Reorder the given list so that elements with tag numbers < the first+-- argument have an index equal to their tag number (if possible)+reorderPrimaryBlocks :: Int -> [(Int, a)] -> [(Int, a)]+reorderPrimaryBlocks n blocks = inner (sortTagged pBlocks) sBlocks [] where+ (pBlocks, sBlocks) = partition (\(tag, _) -> tag < n) blocks+ inner [] sBlocks acc = acc ++ sBlocks+ inner pBlocks [] acc = acc ++ pBlocks+ inner pBlocks@((tag, a) : ps) sBlocks@(s : ss) acc =+ if length acc == tag+ then inner ps sBlocks (acc ++ [(tag, a)])+ else inner pBlocks ss (acc ++ [s])++-- | Recover the primary blocks from a list of @k@ blocks. Each block must be+-- tagged with its number (see the module comments about block numbering)+decode :: FECParams+ -> [(Int, B.ByteString)] -- ^ a list of @k@ blocks and their index+ -> [B.ByteString] -- ^ a list the @k@ primary blocks+decode (FECParams params k n) inblocks+ | length (nub $ map fst inblocks) /= length (inblocks) = error "Duplicate input blocks in FEC decode"+ | any (\f -> f < 0 || f >= n) $ map fst inblocks = error "Invalid block numbers in FEC decode"+ | length inblocks /= k = error "Wrong number of blocks to FEC decode"+ | not (allByteStringsSameLength $ map snd inblocks) = error "Not all inputs to FEC decode are same length"+ | otherwise = unsafePerformIO (do+ let sz = B.length $ snd $ head inblocks+ inblocks' = reorderPrimaryBlocks k inblocks+ presentBlocks = map fst inblocks'+ withForeignPtr params (\cfec -> do+ byteStringsToArray (map snd inblocks') (\src -> do+ b <- createByteStringArray (n - k) sz (\out -> do+ uintCArray presentBlocks (\block_nums -> do+ _decode cfec src out block_nums $ fromIntegral sz))+ let blocks = [0..(n - 1)] \\ presentBlocks+ tagged = zip blocks b+ allBlocks = sortTagged $ tagged ++ inblocks'+ return $ take k $ map snd allBlocks)))++-- | Break a ByteString into @n@ parts, equal in length to the original, such+-- that all @n@ are required to reconstruct the original, but having less+-- than @n@ parts reveals no information about the orginal.+--+-- This code works in IO monad because it needs a source of random bytes,+-- which it gets from /dev/urandom. If this file doesn't exist an+-- exception results+--+-- Not terribly fast - probably best to do it with short inputs (e.g. an+-- encryption key)+secureDivide :: Int -- ^ the number of parts requested+ -> B.ByteString -- ^ the data to be split+ -> IO [B.ByteString]+secureDivide n input+ | n < 0 = error "secureDivide called with negative number of parts"+ | otherwise = withFile "/dev/urandom" ReadMode (\handle -> do+ let inner 1 bs = return [bs]+ inner n bs = do+ mask <- B.hGet handle (B.length bs)+ let masked = B.pack $ B.zipWith xor bs mask+ rest <- inner (n - 1) masked+ return (mask : rest)+ inner n input)++-- | Reverse the operation of secureDivide. The order of the inputs doesn't+-- matter, but they must all be the same length+secureCombine :: [B.ByteString] -> B.ByteString+secureCombine [] = error "Passed empty list of inputs to secureCombine"+secureCombine [a] = a+secureCombine [a, b] = B.pack $ B.zipWith xor a b+secureCombine (a : rest) = B.pack $ B.zipWith xor a $ secureCombine rest++-- | A utility function which takes an arbitary input and FEC encodes it into a+-- number of blocks. The order the resulting blocks doesn't matter so long+-- as you have enough to present to @deFEC@.+enFEC :: Int -- ^ the number of blocks required to reconstruct+ -> Int -- ^ the total number of blocks+ -> B.ByteString -- ^ the data to divide+ -> [B.ByteString] -- ^ the resulting blocks+enFEC k n input = taggedPrimaryBlocks ++ taggedSecondaryBlocks where+ taggedPrimaryBlocks = map (uncurry B.cons) $ zip [0..] primaryBlocks+ taggedSecondaryBlocks = map (uncurry B.cons) $ zip [(fromIntegral k)..] secondaryBlocks+ remainder = B.length input `mod` k+ paddingLength = if remainder >= 1 then (k - remainder) else k+ paddingBytes = (B.replicate (paddingLength - 1) 0) `B.append` (B.singleton $ fromIntegral paddingLength)+ divide a bs+ | B.null bs = []+ | otherwise = (B.take a bs) : (divide a $ B.drop a bs)+ input' = input `B.append` paddingBytes+ blockSize = B.length input' `div` k+ primaryBlocks = divide blockSize input'+ secondaryBlocks = encode params primaryBlocks+ params = fec k n++-- | Reverses the operation of @enFEC@.+deFEC :: Int -- ^ the number of blocks required (matches call to @enFEC@)+ -> Int -- ^ the total number of blocks (matches call to @enFEC@)+ -> [B.ByteString] -- ^ a list of k, or more, blocks from @enFEC@+ -> B.ByteString+deFEC k n inputs+ | length inputs < k = error "Too few inputs to deFEC"+ | otherwise = B.take (B.length fecOutput - paddingLength) fecOutput where+ paddingLength = fromIntegral $ B.last fecOutput+ inputs' = take k inputs+ taggedInputs = map (\bs -> (fromIntegral $ B.head bs, B.tail bs)) inputs'+ fecOutput = B.concat $ decode params taggedInputs+ params = fec k n
+ zfec/fec.c view
@@ -0,0 +1,570 @@+/**+ * zfec -- fast forward error correction library with Python interface+ */++#include "fec.h"++#include <stdio.h>+#include <stdlib.h>+#include <string.h>+#include <assert.h>++/*+ * Primitive polynomials - see Lin & Costello, Appendix A,+ * and Lee & Messerschmitt, p. 453.+ */+static const char*const Pp="101110001";+++/*+ * To speed up computations, we have tables for logarithm, exponent and+ * inverse of a number. We use a table for multiplication as well (it takes+ * 64K, no big deal even on a PDA, especially because it can be+ * pre-initialized an put into a ROM!), otherwhise we use a table of+ * logarithms. In any case the macro gf_mul(x,y) takes care of+ * multiplications.+ */++static gf gf_exp[510]; /* index->poly form conversion table */+static int gf_log[256]; /* Poly->index form conversion table */+static gf inverse[256]; /* inverse of field elem. */+ /* inv[\alpha**i]=\alpha**(GF_SIZE-i-1) */++/*+ * modnn(x) computes x % GF_SIZE, where GF_SIZE is 2**GF_BITS - 1,+ * without a slow divide.+ */+static gf+modnn(int x) {+ while (x >= 255) {+ x -= 255;+ x = (x >> 8) + (x & 255);+ }+ return x;+}++#define SWAP(a,b,t) {t tmp; tmp=a; a=b; b=tmp;}++/*+ * gf_mul(x,y) multiplies two numbers. It is much faster to use a+ * multiplication table.+ *+ * USE_GF_MULC, GF_MULC0(c) and GF_ADDMULC(x) can be used when multiplying+ * many numbers by the same constant. In this case the first call sets the+ * constant, and others perform the multiplications. A value related to the+ * multiplication is held in a local variable declared with USE_GF_MULC . See+ * usage in _addmul1().+ */+static gf gf_mul_table[256][256];++#define gf_mul(x,y) gf_mul_table[x][y]++#define USE_GF_MULC register gf * __gf_mulc_++#define GF_MULC0(c) __gf_mulc_ = gf_mul_table[c]+#define GF_ADDMULC(dst, x) dst ^= __gf_mulc_[x]++/*+ * Generate GF(2**m) from the irreducible polynomial p(X) in p[0]..p[m]+ * Lookup tables:+ * index->polynomial form gf_exp[] contains j= \alpha^i;+ * polynomial form -> index form gf_log[ j = \alpha^i ] = i+ * \alpha=x is the primitive element of GF(2^m)+ *+ * For efficiency, gf_exp[] has size 2*GF_SIZE, so that a simple+ * multiplication of two numbers can be resolved without calling modnn+ */+static void+_init_mul_table(void) {+ int i, j;+ for (i = 0; i < 256; i++)+ for (j = 0; j < 256; j++)+ gf_mul_table[i][j] = gf_exp[modnn (gf_log[i] + gf_log[j])];++ for (j = 0; j < 256; j++)+ gf_mul_table[0][j] = gf_mul_table[j][0] = 0;+}++#define NEW_GF_MATRIX(rows, cols) \+ (gf*)malloc(rows * cols)++/*+ * initialize the data structures used for computations in GF.+ */+static void+generate_gf (void) {+ int i;+ gf mask;++ mask = 1; /* x ** 0 = 1 */+ gf_exp[8] = 0; /* will be updated at the end of the 1st loop */+ /*+ * first, generate the (polynomial representation of) powers of \alpha,+ * which are stored in gf_exp[i] = \alpha ** i .+ * At the same time build gf_log[gf_exp[i]] = i .+ * The first 8 powers are simply bits shifted to the left.+ */+ for (i = 0; i < 8; i++, mask <<= 1) {+ gf_exp[i] = mask;+ gf_log[gf_exp[i]] = i;+ /*+ * If Pp[i] == 1 then \alpha ** i occurs in poly-repr+ * gf_exp[8] = \alpha ** 8+ */+ if (Pp[i] == '1')+ gf_exp[8] ^= mask;+ }+ /*+ * now gf_exp[8] = \alpha ** 8 is complete, so can also+ * compute its inverse.+ */+ gf_log[gf_exp[8]] = 8;+ /*+ * Poly-repr of \alpha ** (i+1) is given by poly-repr of+ * \alpha ** i shifted left one-bit and accounting for any+ * \alpha ** 8 term that may occur when poly-repr of+ * \alpha ** i is shifted.+ */+ mask = 1 << 7;+ for (i = 9; i < 255; i++) {+ if (gf_exp[i - 1] >= mask)+ gf_exp[i] = gf_exp[8] ^ ((gf_exp[i - 1] ^ mask) << 1);+ else+ gf_exp[i] = gf_exp[i - 1] << 1;+ gf_log[gf_exp[i]] = i;+ }+ /*+ * log(0) is not defined, so use a special value+ */+ gf_log[0] = 255;+ /* set the extended gf_exp values for fast multiply */+ for (i = 0; i < 255; i++)+ gf_exp[i + 255] = gf_exp[i];++ /*+ * again special cases. 0 has no inverse. This used to+ * be initialized to 255, but it should make no difference+ * since noone is supposed to read from here.+ */+ inverse[0] = 0;+ inverse[1] = 1;+ for (i = 2; i <= 255; i++)+ inverse[i] = gf_exp[255 - gf_log[i]];+}++/*+ * Various linear algebra operations that i use often.+ */++/*+ * addmul() computes dst[] = dst[] + c * src[]+ * This is used often, so better optimize it! Currently the loop is+ * unrolled 16 times, a good value for 486 and pentium-class machines.+ * The case c=0 is also optimized, whereas c=1 is not. These+ * calls are unfrequent in my typical apps so I did not bother.+ */+#define addmul(dst, src, c, sz) \+ if (c != 0) _addmul1(dst, src, c, sz)++#define UNROLL 16 /* 1, 4, 8, 16 */+static void+_addmul1(register gf*restrict dst, const register gf*restrict src, gf c, size_t sz) {+ USE_GF_MULC;+ const gf* lim = &dst[sz - UNROLL + 1];++ GF_MULC0 (c);++#if (UNROLL > 1) /* unrolling by 8/16 is quite effective on the pentium */+ for (; dst < lim; dst += UNROLL, src += UNROLL) {+ GF_ADDMULC (dst[0], src[0]);+ GF_ADDMULC (dst[1], src[1]);+ GF_ADDMULC (dst[2], src[2]);+ GF_ADDMULC (dst[3], src[3]);+#if (UNROLL > 4)+ GF_ADDMULC (dst[4], src[4]);+ GF_ADDMULC (dst[5], src[5]);+ GF_ADDMULC (dst[6], src[6]);+ GF_ADDMULC (dst[7], src[7]);+#endif+#if (UNROLL > 8)+ GF_ADDMULC (dst[8], src[8]);+ GF_ADDMULC (dst[9], src[9]);+ GF_ADDMULC (dst[10], src[10]);+ GF_ADDMULC (dst[11], src[11]);+ GF_ADDMULC (dst[12], src[12]);+ GF_ADDMULC (dst[13], src[13]);+ GF_ADDMULC (dst[14], src[14]);+ GF_ADDMULC (dst[15], src[15]);+#endif+ }+#endif+ lim += UNROLL - 1;+ for (; dst < lim; dst++, src++) /* final components */+ GF_ADDMULC (*dst, *src);+}++/*+ * computes C = AB where A is n*k, B is k*m, C is n*m+ */+static void+_matmul(gf * a, gf * b, gf * c, unsigned n, unsigned k, unsigned m) {+ unsigned row, col, i;++ for (row = 0; row < n; row++) {+ for (col = 0; col < m; col++) {+ gf *pa = &a[row * k];+ gf *pb = &b[col];+ gf acc = 0;+ for (i = 0; i < k; i++, pa++, pb += m)+ acc ^= gf_mul (*pa, *pb);+ c[row * m + col] = acc;+ }+ }+}++/*+ * _invert_mat() takes a matrix and produces its inverse+ * k is the size of the matrix.+ * (Gauss-Jordan, adapted from Numerical Recipes in C)+ * Return non-zero if singular.+ */+static void+_invert_mat(gf* src, unsigned k) {+ gf c, *p;+ unsigned irow = 0;+ unsigned icol = 0;+ unsigned row, col, i, ix;++ unsigned* indxc = (unsigned*) malloc (k * sizeof(unsigned));+ unsigned* indxr = (unsigned*) malloc (k * sizeof(unsigned));+ unsigned* ipiv = (unsigned*) malloc (k * sizeof(unsigned));+ gf *id_row = NEW_GF_MATRIX (1, k);++ memset (id_row, '\0', k * sizeof (gf));+ /*+ * ipiv marks elements already used as pivots.+ */+ for (i = 0; i < k; i++)+ ipiv[i] = 0;++ for (col = 0; col < k; col++) {+ gf *pivot_row;+ /*+ * Zeroing column 'col', look for a non-zero element.+ * First try on the diagonal, if it fails, look elsewhere.+ */+ if (ipiv[col] != 1 && src[col * k + col] != 0) {+ irow = col;+ icol = col;+ goto found_piv;+ }+ for (row = 0; row < k; row++) {+ if (ipiv[row] != 1) {+ for (ix = 0; ix < k; ix++) {+ if (ipiv[ix] == 0) {+ if (src[row * k + ix] != 0) {+ irow = row;+ icol = ix;+ goto found_piv;+ }+ } else+ assert (ipiv[ix] <= 1);+ }+ }+ }+ found_piv:+ ++(ipiv[icol]);+ /*+ * swap rows irow and icol, so afterwards the diagonal+ * element will be correct. Rarely done, not worth+ * optimizing.+ */+ if (irow != icol)+ for (ix = 0; ix < k; ix++)+ SWAP (src[irow * k + ix], src[icol * k + ix], gf);+ indxr[col] = irow;+ indxc[col] = icol;+ pivot_row = &src[icol * k];+ c = pivot_row[icol];+ assert (c != 0);+ if (c != 1) { /* otherwhise this is a NOP */+ /*+ * this is done often , but optimizing is not so+ * fruitful, at least in the obvious ways (unrolling)+ */+ c = inverse[c];+ pivot_row[icol] = 1;+ for (ix = 0; ix < k; ix++)+ pivot_row[ix] = gf_mul (c, pivot_row[ix]);+ }+ /*+ * from all rows, remove multiples of the selected row+ * to zero the relevant entry (in fact, the entry is not zero+ * because we know it must be zero).+ * (Here, if we know that the pivot_row is the identity,+ * we can optimize the addmul).+ */+ id_row[icol] = 1;+ if (memcmp (pivot_row, id_row, k * sizeof (gf)) != 0) {+ for (p = src, ix = 0; ix < k; ix++, p += k) {+ if (ix != icol) {+ c = p[icol];+ p[icol] = 0;+ addmul (p, pivot_row, c, k);+ }+ }+ }+ id_row[icol] = 0;+ } /* done all columns */+ for (col = k; col > 0; col--)+ if (indxr[col-1] != indxc[col-1])+ for (row = 0; row < k; row++)+ SWAP (src[row * k + indxr[col-1]], src[row * k + indxc[col-1]], gf);+}++/*+ * fast code for inverting a vandermonde matrix.+ *+ * NOTE: It assumes that the matrix is not singular and _IS_ a vandermonde+ * matrix. Only uses the second column of the matrix, containing the p_i's.+ *+ * Algorithm borrowed from "Numerical recipes in C" -- sec.2.8, but largely+ * revised for my purposes.+ * p = coefficients of the matrix (p_i)+ * q = values of the polynomial (known)+ */+void+_invert_vdm (gf* src, unsigned k) {+ unsigned i, j, row, col;+ gf *b, *c, *p;+ gf t, xx;++ if (k == 1) /* degenerate case, matrix must be p^0 = 1 */+ return;+ /*+ * c holds the coefficient of P(x) = Prod (x - p_i), i=0..k-1+ * b holds the coefficient for the matrix inversion+ */+ c = NEW_GF_MATRIX (1, k);+ b = NEW_GF_MATRIX (1, k);++ p = NEW_GF_MATRIX (1, k);++ for (j = 1, i = 0; i < k; i++, j += k) {+ c[i] = 0;+ p[i] = src[j]; /* p[i] */+ }+ /*+ * construct coeffs. recursively. We know c[k] = 1 (implicit)+ * and start P_0 = x - p_0, then at each stage multiply by+ * x - p_i generating P_i = x P_{i-1} - p_i P_{i-1}+ * After k steps we are done.+ */+ c[k - 1] = p[0]; /* really -p(0), but x = -x in GF(2^m) */+ for (i = 1; i < k; i++) {+ gf p_i = p[i]; /* see above comment */+ for (j = k - 1 - (i - 1); j < k - 1; j++)+ c[j] ^= gf_mul (p_i, c[j + 1]);+ c[k - 1] ^= p_i;+ }++ for (row = 0; row < k; row++) {+ /*+ * synthetic division etc.+ */+ xx = p[row];+ t = 1;+ b[k - 1] = 1; /* this is in fact c[k] */+ for (i = k - 1; i > 0; i--) {+ b[i-1] = c[i] ^ gf_mul (xx, b[i]);+ t = gf_mul (xx, t) ^ b[i-1];+ }+ for (col = 0; col < k; col++)+ src[col * k + row] = gf_mul (inverse[t], b[col]);+ }+ free (c);+ free (b);+ free (p);+ return;+}++static int fec_initialized = 0;+static void+init_fec (void) {+ generate_gf();+ _init_mul_table();+ fec_initialized = 1;+}++/*+ * This section contains the proper FEC encoding/decoding routines.+ * The encoding matrix is computed starting with a Vandermonde matrix,+ * and then transforming it into a systematic matrix.+ */++#define FEC_MAGIC 0xFECC0DEC++void+fec_free (fec_t *p) {+ assert (p != NULL && p->magic == (((FEC_MAGIC ^ p->k) ^ p->n) ^ (unsigned long) (p->enc_matrix)));+ free (p->enc_matrix);+ free (p);+}++fec_t *+fec_new(unsigned k, unsigned n) {+ unsigned row, col;+ gf *p, *tmp_m;++ fec_t *retval;++ if (fec_initialized == 0)+ init_fec ();++ retval = (fec_t *) malloc (sizeof (fec_t));+ retval->k = k;+ retval->n = n;+ retval->enc_matrix = NEW_GF_MATRIX (n, k);+ retval->magic = ((FEC_MAGIC ^ k) ^ n) ^ (unsigned long) (retval->enc_matrix);+ tmp_m = NEW_GF_MATRIX (n, k);+ /*+ * fill the matrix with powers of field elements, starting from 0.+ * The first row is special, cannot be computed with exp. table.+ */+ tmp_m[0] = 1;+ for (col = 1; col < k; col++)+ tmp_m[col] = 0;+ for (p = tmp_m + k, row = 0; row < n - 1; row++, p += k)+ for (col = 0; col < k; col++)+ p[col] = gf_exp[modnn (row * col)];++ /*+ * quick code to build systematic matrix: invert the top+ * k*k vandermonde matrix, multiply right the bottom n-k rows+ * by the inverse, and construct the identity matrix at the top.+ */+ _invert_vdm (tmp_m, k); /* much faster than _invert_mat */+ _matmul(tmp_m + k * k, tmp_m, retval->enc_matrix + k * k, n - k, k, k);+ /*+ * the upper matrix is I so do not bother with a slow multiply+ */+ memset (retval->enc_matrix, '\0', k * k * sizeof (gf));+ for (p = retval->enc_matrix, col = 0; col < k; col++, p += k + 1)+ *p = 1;+ free (tmp_m);++ return retval;+}++/* To make sure that we stay within cache in the inner loops of fec_encode()+ and fec_decode(). */+#define STRIDE 8192++void+fec_encode(const fec_t* code, const gf*restrict const*restrict const src, gf*restrict const*restrict const fecs, const unsigned*restrict const block_nums, size_t num_block_nums, size_t sz) {+ unsigned char i, j;+ size_t k;+ unsigned fecnum;+ const gf* p;++ for (k = 0; k < sz; k += STRIDE) {+ size_t stride = ((sz-k) < STRIDE)?(sz-k):STRIDE;+ for (i=0; i<num_block_nums; i++) {+ fecnum=block_nums[i];+ assert (fecnum >= code->k);+ memset(fecs[i]+k, 0, stride);+ p = &(code->enc_matrix[fecnum * code->k]);+ for (j = 0; j < code->k; j++)+ addmul(fecs[i]+k, src[j]+k, p[j], stride);+ }+ }+}++/**+ * Build decode matrix into some memory space.+ *+ * @param matrix a space allocated for a k by k matrix+ */+void+build_decode_matrix_into_space(const fec_t*restrict const code, const unsigned*const restrict index, const unsigned k, gf*restrict const matrix) {+ unsigned char i;+ gf* p;+ for (i=0, p=matrix; i < k; i++, p += k) {+ if (index[i] < k) {+ memset(p, 0, k);+ p[i] = 1;+ } else {+ memcpy(p, &(code->enc_matrix[index[i] * code->k]), k);+ }+ }+ _invert_mat (matrix, k);+}++void+fec_decode(const fec_t* code, const gf*restrict const*restrict const inpkts, gf*restrict const*restrict const outpkts, const unsigned*restrict const index, size_t sz) {+ gf* m_dec = (gf*)alloca(code->k * code->k);+ unsigned char outix=0;+ unsigned char row=0;+ unsigned char col=0;+ build_decode_matrix_into_space(code, index, code->k, m_dec);++ for (row=0; row<code->k; row++) {+ assert ((index[row] >= code->k) || (index[row] == row)); /* If the block whose number is i is present, then it is required to be in the i'th element. */+ if (index[row] >= code->k) {+ memset(outpkts[outix], 0, sz);+ for (col=0; col < code->k; col++)+ addmul(outpkts[outix], inpkts[col], m_dec[row * code->k + col], sz);+ outix++;+ }+ }+}++/**+ * zfec -- fast forward error correction library with Python interface+ * + * Copyright (C) 2007 Allmydata, Inc.+ * Author: Zooko Wilcox-O'Hearn+ * + * This file is part of zfec.+ * + * See README.txt for licensing information.+ */++/*+ * This work is derived from the "fec" software by Luigi Rizzo, et al., the+ * copyright notice and licence terms of which are included below for reference.+ * fec.c -- forward error correction based on Vandermonde matrices 980624 (C)+ * 1997-98 Luigi Rizzo (luigi@iet.unipi.it)+ *+ * Portions derived from code by Phil Karn (karn@ka9q.ampr.org),+ * Robert Morelos-Zaragoza (robert@spectra.eng.hawaii.edu) and Hari+ * Thirumoorthy (harit@spectra.eng.hawaii.edu), Aug 1995+ *+ * Modifications by Dan Rubenstein (see Modifications.txt for + * their description.+ * Modifications (C) 1998 Dan Rubenstein (drubenst@cs.umass.edu)+ *+ * Redistribution and use in source and binary forms, with or without+ * modification, are permitted provided that the following conditions+ * are met:+ *+ * 1. Redistributions of source code must retain the above copyright+ * notice, this list of conditions and the following disclaimer.+ * 2. Redistributions in binary form must reproduce the above+ * copyright notice, this list of conditions and the following+ * disclaimer in the documentation and/or other materials+ * provided with the distribution.+ *+ * THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND+ * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,+ * THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A+ * PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHORS+ * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY,+ * OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,+ * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,+ * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY+ * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR+ * TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT+ * OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY+ * OF SUCH DAMAGE.+ */
+ zfec/fec.h view
@@ -0,0 +1,109 @@+/**+ * zfec -- fast forward error correction library with Python interface+ */++#include <stddef.h>++typedef unsigned char gf;++typedef struct {+ unsigned long magic;+ unsigned k, n; /* parameters of the code */+ gf* enc_matrix;+} fec_t;++#if defined(_MSC_VER)+// actually, some of the flavors (i.e. Enterprise) do support restrict+//#define restrict __restrict+#define restrict+#endif++/**+ * param k the number of blocks required to reconstruct+ * param m the total number of blocks created+ */+fec_t* fec_new(unsigned k, unsigned m);+void fec_free(fec_t* p);++/**+ * @param inpkts the "primary blocks" i.e. the chunks of the input data+ * @param fecs buffers into which the secondary blocks will be written+ * @param block_nums the numbers of the desired check blocks (the id >= k) which fec_encode() will produce and store into the buffers of the fecs parameter+ * @param num_block_nums the length of the block_nums array+ * @param sz size of a packet in bytes+ */+void fec_encode(const fec_t* code, const gf*restrict const*restrict const src, gf*restrict const*restrict const fecs, const unsigned*restrict const block_nums, size_t num_block_nums, size_t sz);++/**+ * @param inpkts an array of packets (size k); If a primary block, i, is present then it must be at index i. Secondary blocks can appear anywhere.+ * @param outpkts an array of buffers into which the reconstructed output packets will be written (only packets which are not present in the inpkts input will be reconstructed and written to outpkts)+ * @param index an array of the blocknums of the packets in inpkts+ * @param sz size of a packet in bytes+ */+void fec_decode(const fec_t* code, const gf*restrict const*restrict const inpkts, gf*restrict const*restrict const outpkts, const unsigned*restrict const index, size_t sz);++#if defined(_MSC_VER)+#define alloca _alloca+#else+#ifdef __GNUC__+#ifndef alloca+#define alloca(x) __builtin_alloca(x)+#endif+#else+#include <alloca.h>+#endif+#endif++/**+ * zfec -- fast forward error correction library with Python interface+ * + * Copyright (C) 2007-2008 Allmydata, Inc.+ * Author: Zooko Wilcox-O'Hearn+ * + * This file is part of zfec.+ * + * See README.txt for licensing information.+ */++/*+ * Much of this work is derived from the "fec" software by Luigi Rizzo, et + * al., the copyright notice and licence terms of which are included below + * for reference.+ * + * fec.h -- forward error correction based on Vandermonde matrices+ * 980614+ * (C) 1997-98 Luigi Rizzo (luigi@iet.unipi.it)+ *+ * Portions derived from code by Phil Karn (karn@ka9q.ampr.org),+ * Robert Morelos-Zaragoza (robert@spectra.eng.hawaii.edu) and Hari+ * Thirumoorthy (harit@spectra.eng.hawaii.edu), Aug 1995+ *+ * Modifications by Dan Rubenstein (see Modifications.txt for + * their description.+ * Modifications (C) 1998 Dan Rubenstein (drubenst@cs.umass.edu)+ *+ * Redistribution and use in source and binary forms, with or without+ * modification, are permitted provided that the following conditions+ * are met:++ * 1. Redistributions of source code must retain the above copyright+ * notice, this list of conditions and the following disclaimer.+ * 2. Redistributions in binary form must reproduce the above+ * copyright notice, this list of conditions and the following+ * disclaimer in the documentation and/or other materials+ * provided with the distribution.+ *+ * THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND+ * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,+ * THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A+ * PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHORS+ * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY,+ * OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,+ * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,+ * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY+ * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR+ * TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT+ * OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY+ * OF SUCH DAMAGE.+ */+