How to reduce the matrix size for Newton Computation for my problem ?
6 views (last 30 days)
Show older comments
This code is computed from the formula :
Basically I convert the Farrow Coefficients to Newton coefficients then compute the filter transfer function .This code is from dsp.FarrowRateConverter to compute Farrow Lagrange coefficients
%%Newton Converter function
function C = designModifiedFarrowCoefficients(N)
if N==0
C = 1;
return
end
% Nondegenerate form, explicitly derived in
% http://www.acoustics.hut.fi/~vpv/publications/vesan_vaitos/ch3_pt2_lagrange.pdf
T = eye(N+1);
for i = 1:N+1
for j = i:N+1
T(i,j) = nchoosek(j-1, i-1)*(floor(N/2))^(j-i);
end
end
C = flipud(T)/fliplr(vander(0:N));
% Force near-zeros and and near-ones to eliminate roundoff error
[Ns, Ds] = rat(C);
C(Ns == 0) = 0;
C(Ns == Ds) = 1;
% Note: below is a shorter alternative
% exponents = (0:N);
% taps = exponents' - floor(N/2);
% C = flipud(inv(taps.^exponents));
end
The solution for >> designModifiedFarrowCoefficients(3)
ans =
-0.1667 0.5000 -0.5000 0.1667
0.5000 -1.0000 0.5000 0
-0.3333 -0.5000 1.0000 -0.1667
0 1.0000 0 0
Its a 4cross4 matrix
Now for transformation computing Td,Td1,Td2 and Tz
which are given as 
function Td1 = compute_Td1(M)
% Create a matrix of binomial coefficients
[I, J] = ndgrid(1:M, 1:M);
binomials = zeros(M, M);
for i = 1:M
binomials(i, 1:i) = arrayfun(@(ii, jj) nchoosek(ii-1, jj-1), I(i, 1:i), J(i, 1:i));
end
% Create a matrix for powers
powers = ((- (M - 1) / 2) .^ (I - J)) .* (J <= I);
% Element-wise multiplication to get Td1
Td1 = binomials .* powers;
end
function Td2 = compute_Td2(M)
% Initialize Td2 with identity
Td2 = eye(M+1);
for n = 2:M+1
for k = 2:n-1
Td2(n, k) = Td2(n-1, k-1) - (n-2) * Td2(n-1, k);
end
end
Td2 = Td2(2:M+1, 2:M+1);
end
function T_D = compute_TD(M)
% Compute Td1 using compute_Td1 function
Td1 = compute_Td1(M);
% Compute Td2 using compute_Td2 function
Td2 = compute_Td2(M);
% Matrix multiplication (cross product)
T_D = Td1 * Td2;
end
function T_z = calculate_Tz(M)
% Initialize the matrix T_z with zeros
T_z = zeros(M, M);
% Fill the matrix T_z using the given formula
for i = 1:M
for j = 1:i % j goes from 1 to i
T_z(i, j) = nchoosek(i-1, j-1) * (-1)^(j+1);
end
end
end
Then the function for calculating Newton Coefficients :
function Q = computeNewtonCoefficients(M)
% Compute Td (T_D in the description)
T_D = compute_TD(M);
% Compute Tz
T_z = calculate_Tz(M);
% Compute the modified Lagrange coefficient matrix P
P = designModifiedFarrowCoefficients(M);
% Compute Q using the transformation Q = (T_D^(-T)) * P * (T_z^(-1))
% Q = inv(T_D') .* P .* inv(T_z);
O = kron(inv(T_D'),P);
Q = kron(O,inv(T_z));
end
solution > computeNewtonCoefficients(3)
ans =
Columns 1 through 10
-0.1667 0 0 0.5000 0 0 -0.5000 0 0 0.1667
-0.1667 0.1667 0 0.5000 -0.5000 0 -0.5000 0.5000 0 0.1667
-0.1667 0.3333 -0.1667 0.5000 -1.0000 0.5000 -0.5000 1.0000 -0.5000 0.1667
0.5000 0 0 -1.0000 0 0 0.5000 0 0 0
0.5000 -0.5000 0 -1.0000 1.0000 0 0.5000 -0.5000 0 0
0.5000 -1.0000 0.5000 -1.0000 2.0000 -1.0000 0.5000 -1.0000 0.5000 0
-0.3333 0 0 -0.5000 0 0 1.0000 0 0 -0.1667
-0.3333 0.3333 0 -0.5000 0.5000 0 1.0000 -1.0000 0 -0.1667
-0.3333 0.6667 -0.3333 -0.5000 1.0000 -0.5000 1.0000 -2.0000 1.0000 -0.1667
0 0 0 1.0000 0 0 0 0 0 0
0 0 0 1.0000 -1.0000 0 0 0 0 0
0 0 0 1.0000 -2.0000 1.0000 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
Columns 11 through 20
0 0 -0.3333 0 0 1.0000 0 0 -1.0000 0
-0.1667 0 -0.3333 0.3333 0 1.0000 -1.0000 0 -1.0000 1.0000
-0.3333 0.1667 -0.3333 0.6667 -0.3333 1.0000 -2.0000 1.0000 -1.0000 2.0000
0 0 1.0000 0 0 -2.0000 0 0 1.0000 0
0 0 1.0000 -1.0000 0 -2.0000 2.0000 0 1.0000 -1.0000
0 0 1.0000 -2.0000 1.0000 -2.0000 4.0000 -2.0000 1.0000 -2.0000
0 0 -0.6667 0 0 -1.0000 0 0 2.0000 0
0.1667 0 -0.6667 0.6667 0 -1.0000 1.0000 0 2.0000 -2.0000
0.3333 -0.1667 -0.6667 1.3333 -0.6667 -1.0000 2.0000 -1.0000 2.0000 -4.0000
0 0 0 0 0 2.0000 0 0 0 0
0 0 0 0 0 2.0000 -2.0000 0 0 0
0 0 0 0 0 2.0000 -4.0000 2.0000 0 0
0 0 -0.1667 0 0 0.5000 0 0 -0.5000 0
0 0 -0.1667 0.1667 0 0.5000 -0.5000 0 -0.5000 0.5000
0 0 -0.1667 0.3333 -0.1667 0.5000 -1.0000 0.5000 -0.5000 1.0000
0 0 0.5000 0 0 -1.0000 0 0 0.5000 0
0 0 0.5000 -0.5000 0 -1.0000 1.0000 0 0.5000 -0.5000
0 0 0.5000 -1.0000 0.5000 -1.0000 2.0000 -1.0000 0.5000 -1.0000
0 0 -0.3333 0 0 -0.5000 0 0 1.0000 0
0 0 -0.3333 0.3333 0 -0.5000 0.5000 0 1.0000 -1.0000
0 0 -0.3333 0.6667 -0.3333 -0.5000 1.0000 -0.5000 1.0000 -2.0000
0 0 0 0 0 1.0000 0 0 0 0
0 0 0 0 0 1.0000 -1.0000 0 0 0
0 0 0 0 0 1.0000 -2.0000 1.0000 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
Columns 21 through 30
0 0.3333 0 0 -0.8333 0 0 2.5000 0 0
0 0.3333 -0.3333 0 -0.8333 0.8333 0 2.5000 -2.5000 0
-1.0000 0.3333 -0.6667 0.3333 -0.8333 1.6667 -0.8333 2.5000 -5.0000 2.5000
0 0 0 0 2.5000 0 0 -5.0000 0 0
0 0 0 0 2.5000 -2.5000 0 -5.0000 5.0000 0
1.0000 0 0 0 2.5000 -5.0000 2.5000 -5.0000 10.0000 -5.0000
0 -0.3333 0 0 -1.6667 0 0 -2.5000 0 0
0 -0.3333 0.3333 0 -1.6667 1.6667 0 -2.5000 2.5000 0
2.0000 -0.3333 0.6667 -0.3333 -1.6667 3.3333 -1.6667 -2.5000 5.0000 -2.5000
0 0 0 0 0 0 0 5.0000 0 0
0 0 0 0 0 0 0 5.0000 -5.0000 0
0 0 0 0 0 0 0 5.0000 -10.0000 5.0000
0 0.1667 0 0 -0.8333 0 0 2.5000 0 0
0 0.1667 -0.1667 0 -0.8333 0.8333 0 2.5000 -2.5000 0
-0.5000 0.1667 -0.3333 0.1667 -0.8333 1.6667 -0.8333 2.5000 -5.0000 2.5000
0 0 0 0 2.5000 0 0 -5.0000 0 0
0 0 0 0 2.5000 -2.5000 0 -5.0000 5.0000 0
0.5000 0 0 0 2.5000 -5.0000 2.5000 -5.0000 10.0000 -5.0000
0 -0.1667 0 0 -1.6667 0 0 -2.5000 0 0
0 -0.1667 0.1667 0 -1.6667 1.6667 0 -2.5000 2.5000 0
1.0000 -0.1667 0.3333 -0.1667 -1.6667 3.3333 -1.6667 -2.5000 5.0000 -2.5000
0 0 0 0 0 0 0 5.0000 0 0
0 0 0 0 0 0 0 5.0000 -5.0000 0
0 0 0 0 0 0 0 5.0000 -10.0000 5.0000
0 0 0 0 -0.1667 0 0 0.5000 0 0
0 0 0 0 -0.1667 0.1667 0 0.5000 -0.5000 0
0 0 0 0 -0.1667 0.3333 -0.1667 0.5000 -1.0000 0.5000
0 0 0 0 0.5000 0 0 -1.0000 0 0
0 0 0 0 0.5000 -0.5000 0 -1.0000 1.0000 0
0 0 0 0 0.5000 -1.0000 0.5000 -1.0000 2.0000 -1.0000
0 0 0 0 -0.3333 0 0 -0.5000 0 0
0 0 0 0 -0.3333 0.3333 0 -0.5000 0.5000 0
0 0 0 0 -0.3333 0.6667 -0.3333 -0.5000 1.0000 -0.5000
0 0 0 0 0 0 0 1.0000 0 0
0 0 0 0 0 0 0 1.0000 -1.0000 0
0 0 0 0 0 0 0 1.0000 -2.0000 1.0000
Columns 31 through 36
-2.5000 0 0 0.8333 0 0
-2.5000 2.5000 0 0.8333 -0.8333 0
-2.5000 5.0000 -2.5000 0.8333 -1.6667 0.8333
2.5000 0 0 0 0 0
2.5000 -2.5000 0 0 0 0
2.5000 -5.0000 2.5000 0 0 0
5.0000 0 0 -0.8333 0 0
5.0000 -5.0000 0 -0.8333 0.8333 0
5.0000 -10.0000 5.0000 -0.8333 1.6667 -0.8333
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
-2.5000 0 0 0.8333 0 0
-2.5000 2.5000 0 0.8333 -0.8333 0
-2.5000 5.0000 -2.5000 0.8333 -1.6667 0.8333
2.5000 0 0 0 0 0
2.5000 -2.5000 0 0 0 0
2.5000 -5.0000 2.5000 0 0 0
5.0000 0 0 -0.8333 0 0
5.0000 -5.0000 0 -0.8333 0.8333 0
5.0000 -10.0000 5.0000 -0.8333 1.6667 -0.8333
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
-0.5000 0 0 0.1667 0 0
-0.5000 0.5000 0 0.1667 -0.1667 0
-0.5000 1.0000 -0.5000 0.1667 -0.3333 0.1667
0.5000 0 0 0 0 0
0.5000 -0.5000 0 0 0 0
0.5000 -1.0000 0.5000 0 0 0
1.0000 0 0 -0.1667 0 0
1.0000 -1.0000 0 -0.1667 0.1667 0
1.0000 -2.0000 1.0000 -0.1667 0.3333 -0.1667
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
How can I reduce this size as it will cause problem in interpolation to 4*4 or in general less .
0 Comments
Answers (0)
See Also
Categories
Find more on Multirate Signal Processing in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
You can also select a web site from the following list
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom(English)
Asia Pacific
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)