bigsplines
index
/.automount/ioasun3/root/disk/users/kybic/work/python_volumes/bigsplines.py

 Module bigsplines provides high level functions
for spline manipulation functions from
helper module BIGsplines (implemented in C)
 
 
Jan Kybic, August 2000
$Id: bigsplines.html,v 1.1 2002/08/26 14:11:16 jkybic Exp $

 
Modules
                      
BIGsplines
Gnuplot
RandomArray
bigpyramid
bigtools
copy
copy_reg
math
multiarray
pickle
string
sys
time
types

 
Classes
                      
MSplineSignal
MSplineSignalList
SplineSignal

 
class MSplineSignal
            Implements multidimensional spline signal defined between
(0..xmax1)x(0..xmax2)...
Initialization:
                s=MSplineSignal(x=signal,c=coefs,xmax=limit,
                               degree=3,shift=floor(degree/2))
 
When neither the signal x, nor its coefficients c are given,
set it to zero. shift is an integer number of spline knots added
to the left (and normally also on the right),
i.e. x(0) coincides with c(shift).
 
The x,c, and xmax are supposed to have the same dimensionality
dim, shift is a scalar.
 
The signal can be evaluated, reduced and expanded by two.  Always,
Reduce(Expand(signal))==signal, the converse conserves the length of
the signal. The expansion is exact, the reduction is currently done
in the sense of L2. There is always one spline coefficient centered at
zero.
 
Note that any spline is representable, no boundary conditions come into
play, except when calculating c from x and when reducing.
 
  
__init__(self, x=None, c=None, xmax=None, degree=3, shift=None)
addSignal(self, s, factor=1.0)
 Adds to myself an MSplineSignal s. My size is not changed,
the signal s is either extended or shortened if needed.
The signal s is optionally premultiplied by a multiplicative
factor 'factor'
copy(self)
 Returns an independent copy of the signal
evalat(self, pts, mask=None)
 Evaluates the signal at arbitrary points pts.
pts should be of dimensionality dim+1 and its
first dimension should be dim, so that pts(:,i) gives
a coordinate of point i
evalderat(self, pts, d, mask=None)
 Evaluates the derivative of the signal with respect to
the component d. Otherwise as evalat
evalfract(self, k=None, n=None)
 Evaluates the signal at the tensor product of
0,1/k[i],2/k[i],3/k[i],...n[i].
n defaults to xmax. n does not have to be an integer
k and n must be vectors of size self.dim
n should not be bigger than xmax
evalint(self, n=None)
 Evaluates signal at integers (0,1,..,n[0])x(0,1,...,n[1])x...
n defaults to floor(xmax)
expand(self)
 Returns a twice-expanded signal, such that fnew(x/2)=fold(x).
The expansion is therefore exact. Note that the number of
added knots (shift) doubles.
getC(self)
getContrib(self, x)
 Returns an array of the size of self.c which contains
the spline values of the corresponding spline functions at point x
getDim(self)
getShape(self)
 Returns the effective shape of the signal, i.e. what would
be returned by evalint
reduce(self)
 Returns a twice-reduced signal, optimal in the L2 sense.
setC(self, c)
setXmax(self, x)
 Change xmax limit. No coefficients are modified.

 
class MSplineSignalList
            A list of MSplineSignals, suitable for representing vector
functions
 
  
__init__(self)
 A primitive constructor. Makes an empty list
addSignal(self, s, factor=1.0)
 Adds to MS..List
copy(self)
 Returns an independent copy of the signal
evalat(self, pts)
 Applies evalat to all components, returns an array of size
(dim,getShape). pts should be (dim,n) as for MSplineSignal.evalat
evalfract(self, k=None, n=None)
 Applies an evalfract to all components, returns an array
expand(self)
 Applies expand to all components.
first(self)
 Returns the first component
fromArray(self, c, n, xmax, degree=3)
 Makes a MSplineSignalList from
an array c splitting it along the first dimension
getContrib(self, x)
 call getContrib
multby(self, x)
 Multiplies all components by a factor of x
reduce(self)
 Applies reduce to all components.
setXmax(self, x)
 apply setXmax to all components
toArray(self)
 Returns the list in the form of an array containing the
coefficients c of the signals.
zero(self, xmax, n, degree=3)
 Makes a MSplineSignalList containing
n signals, capable of representing (multidimensional) interval
0..xmax

 
class SplineSignal
            Implements 1D spline signal defined between 0..xmax
Initialization:
                s=SplineSignal(x=signal,c=coefs,xmax=limit,
                               degree=3,shift=floor(degree/2))
 
When neither the signal x, nor its coefficients c are given,
set it to zero. shift is an integer number of spline knots added
to the left (and normally also on the right),
i.e. x(0) coincides with c(shift)
 
The signal can be evaluated, reduced and expanded by two.  Always,
Reduce(Expand(signal))==signal, the converse conserves the length of
the signal. The expansion is exact, the reduction is currently done
in the sense of L2. There is always one spline coefficient centered at
zero.
 
Note that any spline is representable, no boundary conditions come into
play, except when calculating c from x and when reducing.
 
WARNING: No longer maintained, use MSplineSignal instead.
 
  
__init__(self, x=None, c=None, xmax=None, degree=3, shift=None)
evalat(self, pts, mask=None)
 Evaluates the signal at arbitrary points pts
evalfract(self, k=1, n=None)
 Evaluates the signal at 0,1/k,2/k,3/k,...n.
n defaults to floor(xmax). n does not have to be an integer
evalint(self, n=None)
 Evaluates signal at integers 0,1,..,n.
n defaults to floor(xmax)
expand(self)
 Returns a twice-expanded signal, such that fnew(x/2)=fold(x).
The expansion is therefore exact. Note that the number of
added knots (shift) doubles.
reduce(self)
 Returns a twice-reduced signal, optimal in the L2 sense.

 
Functions
                      
alltrue = reduce(...)
applyseparable(x, oper, passdim=None)
 Apply operation oper on all dimensions of an 
arrays x. Oper should take one argument - a 1D double array and return
an array of the same size. applyseparable returns a double array
 
If passdim is true, the 'oper' is passed a keyword argument
'ndim' giving the dimension being processed.
applyseparablewithcopy(x, oper, targetlen=None, passdim=None)
 Apply operation oper on all dimensions of an 
arrays x. Oper should take one argument - a 1D double array and return
another 1D array not necessarily of the same size.
applyseparablewithcopy returns a double array.
 
If targetlen is not None, it is supposed to be the list of target
dimensions, the respective dimension will be passed to 'oper'.
 
If passdim is true, the 'oper' is passed a keyword argument
'ndim' giving the dimension being processed.
Currently targetlen and ndim are mutually exclusive.
array(...)
array(sequence, typecode=None, copy=1, savespace=0) will return a new array formed from the given (potentially nested) sequence with type given by typecode.  If no typecode is given, then the type will be determined as the minimum type required to hold the objects in sequence.  If copy is zero and sequence is already an array, a reference will be returned.  If savespace is nonzero, the new array will maintain its precision in operations.
bsplnevalfract(k, degree)
 Evaluates a B-spline of degree 'degree' at points
-support+1/k,..,-1/k,0,1/k,..,support-1/k
Returns a 1D float array 'kern' and the offset of the center point
choose(...)
choose(a, (b1,b2,...))
convrect(x, n, bcond=4)
 Convolves 1D vector x with a centered rectangle (Haar) of
length 2*n+1
correctflag(f, d)
 Inverse the order of the axes to pass to the C code for
splineinterpol and similar routines. d is the number of dimensions.
cross_correlate(...)
cross_correlate(a,v)
cumproduct = accumulate(...)
cumsum = accumulate(...)
firspline(c, degree, bcond=4)
 Given B-spline coefficients (c), calculate the function values (x) 
x(k)=sum  c(i) spline_of_degree_N(k-i),   for i=1..n
 
degree=0 constant, 1 linear, 2 quadratic ...
 
Input is scalar. Returns a flat double array
fromstring(...)
fromstring(string, count=None, typecode='l') returns a new 1d array initialized from the raw binary data in string.  If count is not equal to None, the new array will have count elements, otherwise it's size is determined by the size of string.
fspline(x, degree, bcond=4)
 Given values (x) of a function in points 1,2,...n, 
finds coefficients c, such as for all k=1..n
x(k)=sum  c(i) spline_of_degree_N(k-i),   for i=1..n
 
i.e., finds B-spline interpolation
 
degree=0 constant, 1 linear, 2 quadratic ...
 
References : 
M. Unser, A. Aldroubi and M. Eden, Fast B-spline transforms for
continuous image representation and interpolation, IEEE Trans. Pattern Anal.
Machine Intell., vol. 13, pp. 277-285, March 1991.
M. Unser, A. Aldroubi and M. Eden, B-spline signal processing. Part II :
efficient design and applications, IEEE Trans. Signal Processing, vol. 41,
pp. 834-848, February 1993.
 
Returns a flat double array
mconvrect(x, n, bcond=4)
 Convolves multidimensional x with a (multidimensional) rectangle
of edge size 2*n+1
mfirspline(c, degree, bcond=4)
 Finds function values x given B-spline interpolation coefficients
c. Works for any dimensionality. Returns a double array of values
mfspline(x, degree, bcond=4)
 Finds B-spline interpolation coefficients for multidimensional
arrays x. Works for any dimensionality.
Returns a double array of coefs
product = reduce(...)
repeat(...)
repeat(a, n, axis=0)
reshape(...)
reshape(a, (d1, d2, ..., dn)).  Change the shape of a to be an n-dimensional array with dimensions given by d1...dn.  One dimension is allowed to be None.  This dimension will be set to whatever value will make the size of the new array equal the size of the old one.  Note: the size specified for the new array must be exactly equal to the size of the  old one or an error will occur.  This returns a completely new array with the data of the old one copied.  Use a.shape=(...) for no data copying.
sdcgr(ref, cx, cy, degc, derw, wt, mask=None)
 wrapper for BIGsplines.Sdcgr and SdcgrMask.
It is used in bigregister.WarpingProblem.getEgC for fast
calculation of the criterion and gradient
searchsorted = binarysearch(...)
binarysearch(a,v)
sometrue = reduce(...)
splineinterpol(coefs, points, degree, flag=0, mask=None, bcond=4)
 a wrapper to BIGsplines.SplineInterpol. It interpolates a function
given by B-spline coefficients (coefs) at points (points).
Points can by generated by indices(), coefs by mfspline()
flag is passed on to BIGsplines.SplineInterpol to interpolate also
derivatives.
 
If mask is given, the function is only interpolated for those points
for which mask is nonzero
splineinterpolone(coefind, coefsize, points, degree, flag=0, bdw=0, bcond=4)
 a wrapper to BIGsplines.SplineInterpolOne. It interpolates a function
given by B-spline coefficients matrix of size coefsize
knots with exactly one one at point coefind. It evaluates it
at points (points).
Points can by generated by indices(), coefs by mfspline()
flag is passed on to BIGsplines.SplineInterpol to interpolate also
derivatives.
 
Uses mirror boundary conditions
 
bdw is an internal parameter containing the dimension that has been
dealt with
splineinterpoloneraw(coefind, coefsize, points, degree, flag=0, bcond=4)
 a wrapper to BIGsplines.SplineInterpolOne. It interpolates a function
given by B-spline coefficients matrix of size coefsize
knots with exactly one one at point coefind. It evaluates it
at points (points).
Points can by generated by indices(), coefs by mfspline()
flag is passed on to BIGsplines.SplineInterpol to interpolate also
derivatives.
 
Does not consider mirrored coefficients
splinerange(coefind, coefsize, coordsize, shifts, steps, degree, bcond=4)
 returns two vectors l,h containing the lower and upper boundary 
for the support of the B-spline centered at coefind among coefsize
coefs. The returned indices are with respect to an array
 
indices(coordsize,'d')*steps+shifts
sum = reduce(...)
take(...)
take(a, indices, axis=0).  Selects the elements in indices from array a along the given axis.
test_MSlena()
test_MSplineSignal()
test_SplineSignal()
test_bigsplines()
test_evalfract()
test_fir()
test_getContrib()
test_mfspline()
test_splineinterpolone()
time_splineinterpol()
zeros(...)
zeros((d1,...,dn),typecode,savespace) will return a new array of shape (d1,...,dn) and type typecode (default 'l') with all it's entries initialized to zero.  If savespace (default 0) is nonzero the array will be a spacesaver array.

 
Constants
                      ArrayType = <type 'array'>
BufferType = <type 'buffer'>
BuiltinFunctionType = <type 'builtin_function_or_method'>
BuiltinMethodType = <type 'builtin_function_or_method'>
ClassType = <type 'class'>
CodeType = <type 'code'>
Complex = 'D'
Complex0 = 'F'
Complex16 = 'F'
Complex32 = 'F'
Complex64 = 'D'
Complex8 = 'F'
ComplexType = <type 'complex'>
DictType = <type 'dictionary'>
DictionaryType = <type 'dictionary'>
EllipsisType = <type 'ellipsis'>
FileType = <type 'file'>
Float = 'd'
Float0 = 'f'
Float16 = 'f'
Float32 = 'f'
Float64 = 'd'
Float8 = 'f'
FloatType = <type 'float'>
FrameType = <type 'frame'>
FunctionType = <type 'function'>
InstanceType = <type 'instance'>
Int = 'l'
Int0 = '1'
Int16 = 's'
Int32 = 'i'
Int8 = '1'
IntType = <type 'int'>
LambdaType = <type 'function'>
ListType = <type 'list'>
LittleEndian = 1
LongType = <type 'long int'>
MethodType = <type 'instance method'>
ModuleType = <type 'module'>
NoneType = <type 'None'>
PrecisionError = 'PrecisionError'
PyObject = 'O'
SliceType = <type 'slice'>
StringType = <type 'string'>
TracebackType = <type 'traceback'>
TupleType = <type 'tuple'>
TypeType = <type 'type'>
UnboundMethodType = <type 'instance method'>
UnicodeType = <type 'unicode'>
UnsignedInt8 = 'b'
XRangeType = <type 'xrange'>
arraytype = <type 'array'>
defaultbcond = 4
e = 2.7182818284590451
pi = 3.1415926535897931