P      !" # $ % & ' ( ) * + , - . / 0123456789:;<=>?@ABCDEFGHIJKLMNOPQRSTUVWXYZ[\]^_`abcdefghijkl m!n"o"p"q"r"s"t#u#v#w#x$y$z${$|%}&~''''''''(((((()********************************************++,,,,,,,,,,,,,,,,-.......////000112344444445555556778888888888888889999999:;;;; ; < < < <=>>>?@@@@@@@@@AAAA A!A"A#A$A%B&B'C(C)D*D+D,D-D.D/D0D1D2E3E4F5F6F7F8G9G:H;I<J=J>J?J@JAKBLCLDLELFLGLHLILJLKLLLMLNMONP(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafee^xln (1+x), 0 <= x <= 1cos xsin xatan x, -1 < x < 1cosh xsinh xatanh x PQ PQ(c) Matthew Donadio 2002GPLm.p.donadio@ieee.org experimentalportableSafe ,Evaluate a polynomial using Horner's method. Add two polynomials Subtract two polynomials Scale a polynomial Multiply two polynomialsDivide two polynomials2Modulus of two polynomials (remainder of division)2Raise a polynomial to a non-negative integer power&Polynomial substitution y(n) = x(w(n))Polynomial substitution y(n) = x(w(n)) where the coefficients of x are also polynomials.Polynomial derivativePolynomial integrationConvert roots to a polynomial    (c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafegenerates Chebyshev polynomialsNT_N(x)0(c) 1998 Numeric Quest Inc., All rights reservedGPLm.p.donadio@ieee.org experimentalportableSafe#Root finder using Laguerre's methodepsiloniteration limitthe polynomial the rootsRR(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe*Rader's Algorithm using direct convolution'Rader's Algorithm using FFT convolutionx[n]NX[k]x[n]N FFT functionX[k]SS(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeRadix-2 Decimation in Time FFTx[n]N FFT functionX[k](c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe8Prime Factor Algorithm doing row FFT's then column FFT'sx[n]nrowsncols FFT functionX[k]TUVWXYTUVWXY(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe'Goertzel's algorithm for complex inputs1Power via Goertzel's algorithm for complex inputs $Goertzel's algorithm for real inputs!.Power via Goertzel's algorithm for real inputsx[n]kX[k]x[n]k|X[k]|^2 x[n]kX[k]!x[n]k|X[k]|^2 ! ! ! (c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe" Length 2 FFT# Length 3 FFT$ Length 4 FFT"x[n]X[k]#x[n]X[k]$x[n]X[k]"#$"#$"#$ (c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe%x[n]X[k]Z%%%Z (c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe&8Cooley-Tukey algorithm doing row FFT's then column FFT's'8Cooley-Tukey algorithm doing column FFT's then row FFT's&x[n]nrowsncols FFT functionX[k]'x[n]nrowsncols fft functionX[k][\]^_`&'&'&'[\]^_` (c) Matthew Donadio 2002GPLm.p.donadio@ieee.org experimentalportableSafe(Compute the mean of a list Mean(X) = 1/N sum(i=1..N) x_i)Compute the variance of a list Var(X) = sigma^2 % = 1/N-1 sum(i=1..N) (x_i-mu)^2*(Compute the standard deviation of a list " StdDev(X) = sigma = sqrt (Var(X))+'Compute the average deviation of a list % AvgDev(X) = 1/N sum(i=1..N) |x_i-mu|,Compute the skew of a list - Skew(X) = 1/N sum(i=1..N) ((x_i-mu)/sigma)^3-Compute the kurtosis of a list 5 Kurt(X) = ( 1/N sum(i=1..N) ((x_i-mu)/sigma)^4 ) - 3()*+,-()*+,-()*+,-()*+,- (c) Matthew Donadio 2002GPLm.p.donadio@ieee.org experimentalportableSafe.Compute the median of a list/XCompute the center of the list in a more lazy manner and thus halves memory requirement.././././(c) Matthew Donadio 2002GPLm.p.donadio@ieee.org experimentalportableSafe0000(c) Matthew Donadio 2002GPLm.p.donadio@ieee.org experimentalportableSafe1X1X2t2X1X2t3X1X2t123123123Safe 456789:;<=>? 456789:;<=>? 456789:;<=>? 456789:;<=>?(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe@noise white noise@@@(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeAnoise purple noiseAAA(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeBnoise brown noiseBBB(c) Matt Harden 1999GPLm.p.donadio@ieee.org experimentalportableSafeCabcdefghijklmnopqDECDECDECabcdefghijklmnopqDEa8b8(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe F32 bits in [0,1]G32 bits in [0,1)H32 bits in (0,1]I32 bits in (0,1)J-53 bits in [0,1], ie 64-bit IEEE 754 in [0,1]K-53 bits in [0,1), ie 64-bit IEEE 754 in [0,1)L53 bits in (0,1]M53 bits in (0,1)N!transforms uniform [0,1] to [a,b] FXUGXUHXUIXUJXUKXULXUMXUNabUU' FGHIJKLMN FGHIJKLMN FGHIJKLMN(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeOEGenerates a list of poisson random variables from a list of uniforms.r7Split after every element that satisfies the predicate.!A candidate for a Utility module.O3lambda - expectation value, should be non-negative.4uniformly distributed values from the interval [0,1]Poisson distributed outputsrPQOPQOPQOrPQ(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeRGGenerates a list of geometric random variables from a list of uniformsRpUXRRR(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeSiGenerates a list of gamma random variables from a list of uniforms via the inverse transformation methodSnlambdaUXSSS(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeToGenerates a list of exponential random variables from a list of uniforms via the inverse transformation method  F(x) = 1 - exp(-lambda*x)  F^-1(x) = -log(1 - x) / lambdaTlambdaUXTTT(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeUFGenerates a list of binomial random variables from a list of uniformsUnpUXUUU(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeV*Calculates the Chebyshev approximation to f(x) over [a,b]W)Evaluates the Chebyshev approximation to f(x) over [a,b] at xVf(x)abNc_nWc_nabxf(x)VWVWVW(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeX)Type for results of the simplex algorithm\(The simplex algorithm for standard form: min c'xwhere Ax = b, x >= 0 a!(0,0) = -z a!(0,j) = c' a!(i,0) = ba!(i,j) = A_ij]The two-phase simplex algorithmXYZ[stuvw\stating tableausolutionxyz]stating tableausolution{|XYZ[\]XYZ[\] XYZ[stuvw\xyz]{|(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe`'Matrix-matrix multiplication: A x B = Ca'Matrix-vector multiplication: A x b = cbTranspose of a matrixc5Hermitian transpose (conjugate transpose) of a matrix`ABCaAbcbAA^TcAA^H`abc`abc`abc(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafedlevinson takes an array, r, of autocorrelation values, and a model order, p, and returns an array, a, of the model estimate and rho, the noise power.drp(a,rho)ddd(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafee&LU decomposition via Crout's Algorithmf%Solution to Ax=b via LU decompositiong/Improve a solution to Ax=b via LU decompositionh%Matrix inversion via LU decompositioni,Determinant of a matrix via LU decompositionjLU solver using original matrixk!determinant using original matrixeALU(A)fLU(A)bxgALU(A)bxx'hAA^-1iLU(A)det(A)jAbxkAdet(A)efghijkefghijkefghijk (c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafelAbxlll!(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafemJThis is the simple phase unwrapping algorithm from Oppenheim and Schafer.mepsilonARGargmmm"(c) Matthew Donadio 1998,2003GPLm.p.donadio@ieee.org experimentalportableSafennk creates a sine wave with normalized frequency wn (numerically controlled oscillator, or NCO) using the recurrence relation y[n] = 2cos(wn)*y[n-1] - y[n-2]. Eventually, cumulative errors will creep into the data. This is unavoidable since performing AGC on this type of real data is hard. The good news is that the error is small with floating point data.oo mixes (multiplies) x by a real sine wave with normalized frequency wn. This is usually called an NCOM: Numerically Controlled Oscillator and Modulator.qqa returns an infinite list representing a complex phasor with a phase step of wn radians, ie a quadrature nco with normalized frequency wn radians/sample. Since Haskell uses lazy evaluation, rotate will only be computed once, so this NCO uses only one sin and one cos for the entire list, at the expense of 4 mults, 1 add, and 1 subtract per point.rr mixes the complex input x with a quardatue nco with normalized frequency wn radians/sample using complex multiplies (perform a complex spectral shift)ss| mixes the complex input x with a quadrature nco with normalized frequency wn radians/sample in quadrature (I/Q modulation)nwphiyowphixypqwphiyrwphixy}swphixynopqrsnoqrspnopqr}s#(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafet all zerosusingle impulsev unit stepwramptuvwtuvwtuvwtuvw$(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafexCascade of functions, eg ,cascade [ f1, f2, f3 ] x == (f3 . f2 . f1) xy Gain node y[n] = a * x[n]z Bias node y[n] = x[n] + a{ Adder node z[n] = x[n] + y[n]x f_n(x)x[n]y[n]yax[n]y[n]zax[n]y[n]{x[n]y[n]z[n]xyz{xyz{xyz{%(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe| Implementation of Prony's method|pqg[n](b,a)|||&(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe} Performs the matched-z transform}T_s(b,a)(b',a')}}}'(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafePerforms the bilinear transform!Function for frequency prewarping~T_s(b,a)(b',a')w_cT_sW_c~~~((c) Matthew Donadio 1998GPLm.p.donadio@ieee.org experimentalportableSafeLowpass filterHighpass filterBandpass filterBandstop filterMultiband filterRaised-cosine filter ~wcMh[n]wcMh[n]wlwuMh[n]wlwuMh[n] mags wMh[n]wsbetaMh[n] ~)(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafedesigns smooth FIR filterspqh[n]*(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeN takes the continuous impluse response function (one of the functions below, f-) and number of points in the interpolation, p, time shifts it by xh, samples it, and creates an array with the interpolation coeficients that can be used as a FIR filter.,fpxh[n],,,+(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe8Implements the following function, which is a FIR filter y[n] = sum(k=0,M) h[k]*x[n-k]We implement the fir function with five helper functions, depending on the type of the filter. In the following functions, we use the O&S convention that m is the order of the filter, which is equal to the number of taps minus one.h[n]x[n]y[n],(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe =This is an integrator when a==1, and a leaky integrator when  0 < a < 1. y[n] = a * y[n-1] + x[n]First order section, DF1 v[n] = b0 * x[n] + b1 * x[n-1] y[n] = v[n] - a1 * y[n-1]First order section, DF2 w[n] = -a1 * w[n-1] + x[n] y[n] = b0 * w[n] + b1 * w[n-1]First order section, DF2T v0[n] = b0 * x[n] + v1[n-1]  y[n] = v0[n] v1[n] = -a1 * y[n] + b1 * x[n](Direct Form I for a second order section ,v[n] = b0 * x[n] + b1 * x[n-1] + b2 * x[n-2] 'y[n] = v[n] - a1 * y[n-1] - a2 * y[n-2]2Direct Form II for a second order section (biquad) (w[n] = -a1 * w[n-1] - a2 * w[n-2] + x[n] ,y[n] = b0 * w[n] + b1 * w[n-1] + b2 * w[n-2]4Transposed Direct Form II for a second order section v0[n] = b0 * x[n] + v1[n-1]  y[n] = v0[n] (v1[n] = -a1 * y[n] + b1 * x[n] + v2[n-1] v2[n] = -a2 * y[n] + b2 * x[n]Direct Form I IIR v[n] = sum(k=0..M) b_k*x[n-k] $y[n] = v[n] - sum(k=1..N) a_k*y[n-k]v[n] is calculated with Direct Form II IIR $w[n] = x[n] - sum(k=1..N) a_k*w[n-k] y[n] = sum(k=0..M) b_k*w[n-k]ax[n]y[n]a_1b_0b_1x[n]y[n]a_1b_0b_1x[n]y[n]a_1b_0b_1x[n]y[n]a_1a_2b_0b_1b_2x[n]y[n]a_1a_2b_0b_1b_2x[n]y[n]a_1a_2b_0b_1b_2x[n]y[n](b,a)x[n]y[n](b,a)x[n]y[n]-(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafePolyphase interpolatorLh[n]x[n]y[n].(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeLowpass to lowpass:  s --> s/wcLowpass to highpass:  s --> wc/sLowpass to bandpass:  s --> (s^2 + wl*wu) / (s(wu-wl))Lowpass to bandstop:  s --> (s(wu-wl)) / (s^2 + wl*wu) wc(b,a)(b',a')wc(b,a)(b',a')wlwu(b,a)(b',a')wlwu(b,a)(b',a') /(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeLowpass to lowpass:  z^-1 --> (z^-1 - a)/(1 - a*z^-1)Lowpass to Highpass: !z^-1 --> -(z^-1 + a)/(1 + a*z^-1)Lowpass to Bandpass: z^-1 -->Lowpass to Bandstop: z^-1 -->theta_pomega_p(b,a)(b',a')theta_pomega_p(b,a)(b',a')theta_pomega_p1omega_p2(b,a)(b',a')theta_pomega_p1omega_p2(b,a)(b',a')0(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe&Generates Butterworth filter prototype$Generates Chebyshev filter prototype,Generates Inverse Chebyshev filter prototypeN(b,a)epsilonN(b,a)epsilonN(b,a)1(c) Matthew Donadio 2002GPLm.p.donadio@ieee.org experimentalportableSafeComplex test dataReal test data2(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeThe Quinn-Fernandes algorithmyinitial w estimatew3(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe(Discrete frequency periodigram maximizerX[k]kw4(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSaferaw cross-correllationbiased cross-correllationunbiased cross-correllationraw auto-correllationbiased auto-correllationunbiased auto-correllation xykR_xy[k]xyk R_xy[k] / NxykR_xy[k] / (N-k)xkR_xx[k]xk R_xx[k] / NxkR_xx[k] / (N-k) 5(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSaferaw cross-covariance-We define covariance in terms of correlation._Cxy(X,Y) = E[(X - E[X])(Y - E[Y])] = E[XY] - E[X]E[Y] = Rxy(X,Y) - E[X]E[Y]raw auto-covariance]Cxx(X,X) = E[(X - E[X])(X - E[X])] = E[XX] - E[X]E[X] = Rxy(X,X) - E[X]^2biased cross-covarianceunbiased cross-covariancebiased auto-covarianceunbiased auto-covariancexykC_xy[k]xkC_xx[k]xyk C_xy[k] / NxykC_xy[k] / (N-k)xk C_xx[k] / NxkC_xx[k] / (N-k)6(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeTHIS DOES NOT WORK7(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeconv convolves two finite sequences8(c) Matthew Donadio 1998GPLm.p.donadio@ieee.org experimentalportableSafe  generates a list of values linearly spaced between specified start and end values (array will include both start and end values). 2linspace 0.0 1.0 5 == [ 0.0, 0.25, 0.5, 0.75 1.0 ] generates a list of values logarithmically spaced between the values 10 ** start and 10 ** end (array will include both start and end values). 3logspace 0.0 1.0 4 == [ 1.0, 2.1544, 4.6416, 10.0 ] is the unit delay function, eg, $delay1 [ 1, 2, 3 ] == [ 0, 1, 2, 3 ]$ is the n sample delay function, eg, +delay 3 [ 1, 2, 3 ] == [ 0, 0, 0, 1, 2, 3 ] downsample# throws away every n'th sample, eg, 0downsample 2 [ 1, 2, 3, 4, 5, 6 ] == [ 1, 3, 5 ]upsample+ inserts n-1 zeros between each sample, eg, .upsample 2 [ 1, 2, 3 ] == [ 1, 0, 2, 0, 3, 0 ]upsampleAndHold$ replicates each sample n times, eg, >upsampleAndHold 3 [ 1, 2, 3 ] == [ 1, 1, 1, 2, 2, 2, 3, 3, 3 ]Bmerges elements from two lists into one list in an alternating way ;interleave [0,1,2,3] [10,11,12,13] == [0,10,1,11,2,12,3,13]1split a list into two lists in an alternating way /uninterleave [1,2,3,4,5,6] == ([1,3,5],[2,4,6])It's a special case of GO.%pad a sequence with zeros to length n )pad [ 1, 2, 3 ] 6 == [ 1, 2, 3, 0, 0, 0 ] generates a  if the given condition holds7Computes the square of the Euclidean norm of a 2D pointTPower with fixed exponent type. This eliminates warnings about using default types.89(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeQuinn's First Estimator (FCI1)Quinn's Second Estimator (FCI2)Quinn's Third Estimator (FCI3)Eric Jacobsen's EstimatorMacLeod's Three Point EstimatorMacLeod's Three Point EstimatorRife and Vincent's EstimatorX[k]kwX[k]kwX[k]kwX[k]kwX[k]kwX[k]kwX[k]kw:(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe(Pisarenko's method for a single sinusoidxw;(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe=The weighted linear predictor form of the frequency estimator)WLP using Lank, Reed, and Pollon's windowWLP using kay's window(WLP using Lovell and Williamson's window 2WLP using Clarkson, Kootsookos, and Quinn's windowwindowzwzwzwzw zrhosigmaw   <(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe DComputes an AR(p) model estimate from x using the Yule-Walker method CComputes an AR(p) model estimate from x using the covariance method LComputes an AR(p) model estimate from x using the modified covariance method >Computes an AR(p) model estimate from x using the Burg' method xp(a,rho) xp(a,rho) xp(a,rho) xp(a,rho)            =(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeComputes an MA(q) model estimate from x using the Durbin's method where l is the order of the AR process used in the algorithmxql(a,rho)>(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe$Butterworth filter response function"Chebyshev filter response function*Inverse Chebyshev filter response function Note that w_c. is a property of the stopband for this filterNw_cw |H_c(w)|^2Nepsilonw_cw |H_c(w)|^2Nepsilonw_cw |H_c(w)|^2?(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeFilter shaprening routineh[n]-function that implements the sharpened filter@(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe    A(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe)Generates lowpass Butterworth IIR filters )Generates lowpass Butterworth IIR filters!'Generates lowpass Chebyshev IIR filters"'Generates lowpass Chebyshev IIR filters#/Generates lowpass Inverse Chebyshev IIR filters$/Generates lowpass Inverse Chebyshev IIR filters (wp,dp)(ws,ds)(b,a) (wp,dp)(ws,ds)(b,a)!(wp,dp)(ws,ds)(b,a)"(wp,dp)(ws,ds)(b,a)#(wp,dp)(ws,ds)(b,a)$(wp,dp)(ws,ds)(b,a)  !"#$ !# "$  !"#$B(c) Matthew Donadio 1998GPLm.p.donadio@ieee.org experimentalportableSafe%CIC interpolator&CIC interpolator%RMNx[n]y[n]&RMNx[n]y[n]%&%&%&C(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe'Halfband interpolator(Halfband decimator'h[n]x[n]y[n](h[n]x[n]y[n]'('('(D(c) Matthew Donadio 1998GPLm.p.donadio@ieee.org experimentalportableSafe )Applys a window, w, to a sequence x*rectangular window+Bartlett window,Hanning window-Hamming window.Blackman window/Generalized Hamming window0rectangular window1rectangular window)w[n]x[n] w[n] * x[n]*Mw[n]+Mw[n],Mw[n]-Mw[n].Mw[n]/alphaMw[n]0betaMw[n]1Mw[n] )*+,-./01 )*+,-.0/1)*+,-./01E(c) Matthew Donadio 1998GPLm.p.donadio@ieee.org experimentalportableSafe2Designs a lowpass Kaiser filter3 Designs a highpass Kaiser filter2wpwsdpdsh[n]3wpwsdpdsh[n]232323PSafe )*+,-./01F(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe4YNormal random variables via the Central Limit Theorm (not explicity given, but see Ross)FIf mu=0 and sigma=1, then this will generate numbers in the range [-n2,n2]5MNormal random variables via the Box-Mueller Polar Method (Ross, pp 450--452)If mu=0 and sigma=1, then this will generate numbers in the range [-8.57,8.57] assuing that the uniform RNG is really giving full precision for doubles.6/Acceptance-Rejection Method (Ross, pp 448--450)If mu=0 and sigma=1, then this will generate numbers in the range [-36.74,36.74] assuming that the uniform RNG is really giving full precision for doubles.7>Ratio Method (Kinderman-Monahan) (Knuth, v2, 2ed, pp 125--127)CIf mu=0 and sigma=1, then this will generate numbers in the range  ? -1e15,1e15L assuming that the uniform RNG is really giving full precision for doubles.4Number of uniforms to sum (mu,sigma)UX5 (mu,sigma)UX6 (mu,sigma)UX7 (mu,sigma)UX456745674567G(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe8Kellet's filter9Voss's algorithmMUNTESTED, but the algorithm looks like it is working based on my hand tests.8noise pinked noise9number of octaves to sumnoise pinked noise898989H(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe:#Radix-2 Decimation in Frequency FFT:x[n]N FFT functionX[k]:::I(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe;#Radix-4 Decimation in Frequency FFT;x[n]N FFT functionX[k];;;J(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafe<{This is the driver routine for calculating FFT's. All of the recursion in the various algorithms are defined in terms of <.=;Inverse FFT, including scaling factor, defined in terms of <>hThis is the algorithm for computing 2N-point real FFT with an N-point complex FFT, defined in terms of <?rThis is the algorithm for computing a 2N-point real inverse FFT with an N-point complex FFT, defined in terms of =@[Algorithm for 2 N-point real FFT's computed with N-point complex FFT, defined in terms of <<x[n]X[k]=X[k]x[n]>x[n]X[k]?X[k]x[n]@x1[n]x2[n] (X1[k],X2[k])<=>?@<=>?@<=>?@K(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeA fast_conv7 convolves two finite sequences using DFT relationshipsAAAAQSafeFor small values of narrowj (large width) it is faster to compute the spectrum of the Gaussian and apply the Fourier transform on it. For small c: the convergence is slow, but the result will be close to recip (sqrt c). array sizereciprocal of widthX[k] array size6reciprocal of width, 1 is the width of the eigenvectorX[k]L(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeBCDEFGHIJKLM BCDEFGHIJKLM BCDEFGHIJKLMBCDEFGHIJKLMM(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeN Sliding FFTNNx[n] X[k]NNNN(c) Matthew Donadio 2003GPLm.p.donadio@ieee.org experimentalportableSafeO'Split-Radix Decimation in Frequency FFTOx[n]N FFT functionX[k]OOORSTUVWXYZ[\]^_`abcdefghijklmnopqrs t u v w x y z { | } ~   !""""""####$$$$%&''''''''(((((()******************************************** + +, , , ,,,,,,,,,,,,,-..... .!."/#/$/%/&0'0(0)1*1+2,3-4.4/4041424345455565758596:7;78<8=8>8?8@8A8B8C8D8E8F8G8H8I8J9K9L9M9N9O9P9Q:R;S;T;U;V;W<X<Y<Z<[=\>]>^>_?`@@@a@b@c@d@e@f@gAhAiAjAkAlAmAnAoApBqBrCsCtDuDvDwDxDyDzD{D|D}E~EFFFFGGHIJJJJJKLLLLLLLLLLLLMN       "((((())))++++++++++++++++++++++++++,,,,,,,,,,,-..23334447779:;;A'AABBBCDDDD D D D D DDEEEEEEFFGGOGJJJQQQJQQ Q!Q"Q#Q$Q%Q&Q'Q(Q)L*LL+LL,L-L.M/M01"dsp-0.2.3.1-3uBipKNU9sV5aIe7q3ueIfPolynomial.MaclaurinPolynomial.BasicPolynomial.ChebyshevPolynomial.RootsNumeric.Transform.Fourier.RaderNumeric.Transform.Fourier.R2DITNumeric.Transform.Fourier.PFA"Numeric.Transform.Fourier.Goertzel!Numeric.Transform.Fourier.FFTHardNumeric.Transform.Fourier.DFTNumeric.Transform.Fourier.CTNumeric.Statistics.MomentNumeric.Statistics.MedianNumeric.Statistics.CovarianceNumeric.Statistics.TTestNumeric.Special.TrigonometricNumeric.Random.Spectrum.WhiteNumeric.Random.Spectrum.PurpleNumeric.Random.Spectrum.Brown Numeric.Random.Generator.MT19937#Numeric.Random.Distribution.Uniform#Numeric.Random.Distribution.Poisson%Numeric.Random.Distribution.Geometric!Numeric.Random.Distribution.Gamma'Numeric.Random.Distribution.Exponential$Numeric.Random.Distribution.BinomialNumeric.Approximation.ChebyshevMatrix.Simplex Matrix.MatrixMatrix.Levinson Matrix.LUMatrix.Cholesky DSP.UnwrapDSP.Source.OscillatorDSP.Source.Basic DSP.FlowgraphDSP.Filter.IIR.PronyDSP.Filter.IIR.MatchedzDSP.Filter.IIR.BilinearDSP.Filter.FIR.TapsDSP.Filter.FIR.SmoothDSP.Filter.FIR.PolyInterpDSP.Filter.FIR.FIRDSP.Filter.IIR.IIRDSP.Multirate.PolyphaseDSP.Filter.Analog.TransformDSP.Filter.IIR.TransformDSP.Filter.Analog.PrototypeDSP.Estimation.Spectral.KayData'DSP.Estimation.Frequency.QuinnFernandesDSP.Estimation.Frequency.PerMaxDSP.CorrelationDSP.CovarianceDSP.Estimation.Spectral.ARMADSP.Convolution DSP.BasicDSP.Estimation.Frequency.FCI"DSP.Estimation.Frequency.PisarenkoDSP.Estimation.Frequency.WLPDSP.Estimation.Spectral.ARDSP.Estimation.Spectral.MADSP.Filter.Analog.ResponseDSP.Filter.FIR.SharpenDSP.Filter.IIR.CookbookDSP.Filter.IIR.DesignDSP.Multirate.CICDSP.Multirate.Halfband DSP.WindowDSP.Filter.FIR.Kaiser"Numeric.Random.Distribution.NormalNumeric.Random.Spectrum.PinkNumeric.Transform.Fourier.R2DIFNumeric.Transform.Fourier.R4DIFNumeric.Transform.Fourier.FFTDSP.FastConvolution"Numeric.Transform.Fourier.FFTUtils$Numeric.Transform.Fourier.SlidingFFTNumeric.Transform.Fourier.SRDIFsplitDSP.Filter.FIR.Window%Numeric.Transform.Fourier.Eigensystempolyexppolyln1polycospolysinpolyatanpolycoshpolysinh polyatanhpolyevalpolyadd polyAddScalarpolysub polyscalepolymult polymultAltpolydivpolymodpolypow polysubst polysubstAlt polyPolySubst polyderiv polyinteg roots2polychebyroots fft_rader1 fft_rader2 fft_r2ditfft_pfa cgoertzelcgoertzel_power rgoertzelrgoertzel_powerfft'2fft'3fft'4dftfft_ct1fft_ct2meanvarstddevavgdevskewkurtosismedian medianFastcovttesttutesttptestcscseccotacscasecacotcschsechcothacschasechacothwhitepurplebrownWgenrandtest uniform32cc uniform32co uniform32oc uniform32oo uniform53cc uniform53co uniform53oc uniform53oouniformpoissontestHead geometricgammaexponential_invbinomial cheby_approx cheby_evalSimplex Unbounded InfeasibleOptimalsimplextwophase $fReadSimplex $fShowSimplexmm_multmv_multm_transm_hermitlevinsonlulu_solveimproveinverselu_detsolvedetcholeskyunwrapnconcomagcquadrature_nco complex_ncomquadrature_ncomzerosimpulsesteprampcascadegainbiasadderpronymatchedzzmzpstep1step2step3step4bilinearprewarplpfhpfbpfbsfmbfrc smoothfirmkcoef bspline_1p0o bspline_2p1o bspline_4p3o bspline_6p5o lagrange_4p3o lagrange_6p5o hermite_4p3o hermite_6p3o hermite_6p5o sndosc_4p5o sndosc_6p5o watte_4p2oparabolic2x_4p2ooptimal_2p3o2xoptimal_2p3o4xoptimal_2p3o8xoptimal_2p3o16xoptimal_2p3o32xoptimal_4p2o2xoptimal_4p2o4xoptimal_4p2o8xoptimal_4p2o16xoptimal_4p2o32xoptimal_4p3o2xoptimal_4p3o4xoptimal_4p3o8xoptimal_4p3o16xoptimal_4p3o32xoptimal_4p4o2xoptimal_4p4o4xoptimal_4p4o8xoptimal_4p4o16xoptimal_4p4o32xoptimal_6p4o2xoptimal_6p4o4xoptimal_6p4o8xoptimal_6p4o16xoptimal_6p4o32xoptimal_6p5o2xoptimal_6p5o4xoptimal_6p5o8xoptimal_6p5o16xoptimal_6p5o32xfir integratorfos_df1fos_df2fos_df2t biquad_df1 biquad_df2 biquad_df2tiir_df1iir_df2xtytf1f2f3f4f5 poly_interpa_lp2lpa_lp2hpa_lp2bpa_lp2bs substitutepropSubstituteRecippropSubstituteAltd_lp2lpd_lp2hpd_lp2bpd_lp2bs butterworth chebyshev1 chebyshev2xcxrqfpermaxrxyrxy_brxy_urxxrxx_brxx_ucxycxxcxy_bcxy_ucxx_bcxx_u arma_myweconvlinspacelogspacedelay1delay downsample downsampleRecupsample upsampleRecupsampleAndHold interleave uninterleavepadtoMaybenorm2sqr^!quinn1quinn2quinn3jacobsenmacleod3macleod5rv pisarenkowlplrpkaylwckqar_ywar_covar_mcovar_burg ma_durbin butterworth_H chebyshev1_H chebyshev2_Hsharpenbpf_csgbpf_cpgnotchapf peakingEQlowShelf highShelfpoly2iirbutterworthLowpassbutterworthHighpassbutterworthBandpass mkButterworthchebyshev1Lowpass mkChebyshev1chebyshev2Lowpass mkChebyshev2cic_interpolate cic_decimate hb_interphb_decimwindow rectangularbartletthanninghammingblackman gen_hammingkaiserparzen kaiser_lpf kaiser_hpf normal_clt normal_bm normal_arnormal_rkelletvoss fft_r2dif fft_r4diffftifftrfftirfftr2fft fast_convfft_magfft_db fft_phasefft_grdfft_inforfft_magrfft_db rfft_phaserfft_grd rfft_infowrite_fft_infowrite_rfft_infosfft fft_srdififacsinverseslaguerre generator pfa_index_map find_inverserowscols flatten_rows flatten_colsdft' ct_index_map1 ct_index_map2.<<..>>.parmNparmMparmA upperMask lowerMasktemperingMaskBtemperingMaskCtemperingShiftUtemperingShiftStemperingShiftTtemperingShiftLsgenrandmag01 temperingrand segmentAfterepspivotchooseqchooseprefinepaddartcolsumdelartgettabcostquadrature_multindexeslpf_taphpf_tapmbf_taprc_tap normalizeexpandexpand'reflectfir0fir'0fir'1fir'2fir'3fir'4 isFIRType1 isFIRType2 isFIRType3 isFIRType4htyt'h1h2h3h4y1y2y3y4y1'y2'y3'y4' integrator'fos_df1'fos_df2' fos_df2t'df1df2df2tiir'df1iir'df2atbtmkpoly substituteAltqf'calc_xcalc_ypermax'rtbaseGHC.BaseJustlog10rsskayWin processArraybutterworthParamswcapply integratecomb mkhalfband bartlett'hanning'hamming' blackman' gen_hamming'kaiser'i0i0'parzen' makeArraycalc_wccalc_dwcalc_A calc_betacalc_M ceilingEvenadjust normalDistadd mkOctaveschoose1choose2 choose_factorgaussianFourierwrappedGaussian gaussianPlainexample0example1exampleAlgebraicDoubleexampleAlgebraicQuadroexampleAlgebraicHalfexampleAlgebraicFourthexampleAlgebraic2exampleAlgebraic3exampleAlgebraic4exampleAlgebraic5magsqdot hPrintIndex write_cvector write_rvectorsfft'enough