MATLAB: How to convert a code from C language to Matlab

code conversionMATLAB

#include <math.h>
#include <stdio.h>
/*********************************************************************
*
*
* Name: mag_est.c
*
* Synopsis:
*
*   Demonstrates and tests the "Alpha * Min + Beta * Max" magnitude
*   estimation algorithm.
*
* Description:
*
*   This program demonstrates the "Alpha, Beta" algorithm for 
*   estimating the magnitude of a complex number.  Compared to
*   calculating the magnitude directly using sqrt(I^2 + Q^2), this
*   estimation is very quick.
*
*   Various values of Alpha and Beta can be used to trade among RMS
*   error, peak error, and coefficient complexity.  This program
*   includes a table of the most useful values, and it prints out the
*   resulting RMS and peak errors.
*
* Copyright 1999  Grant R. Griffin
*
*                    The Wide Open License (WOL)
*
* Permission to use, copy, modify, distribute and sell this software
* and its documentation for any purpose is hereby granted without
* fee, provided that the above copyright notice and this license
* appear in all source copies. THIS SOFTWARE IS PROVIDED "AS IS"
* WITHOUT EXPRESS OR IMPLIED WARRANTY OF ANY KIND. See 
* http://scopeplot.com/dg/wol.htm for more information.
*
*********************************************************************/
/********************************************************************/
double alpha_beta_mag(double alpha, double beta, double inphase,
                      double quadrature)
{
   /* magnitude ~= alpha * max(|I|, |Q|) + beta * min(|I|, |Q|) */
   double abs_inphase = fabs(inphase);
   double abs_quadrature = fabs(quadrature);
   if (abs_inphase > abs_quadrature) {
      return alpha * abs_inphase + beta * abs_quadrature;
   } else {
      return alpha * abs_quadrature + beta * abs_inphase;
   }
}
/*********************************************************************/
double decibels(double linear)
{
   #define SMALL 1e-20
   if (linear <= SMALL) {
      linear = SMALL;
   }
   return 20.0 * log10(linear);
}
/*********************************************************************/
void test_alpha_beta(char *name, double alpha, double beta,
                     int num_points)
{
   #define PI 3.141592653589793
   int ii;
   double phase, real, imag, err, avg_err, rms_err;
   double peak_err = 0.0;
   double sum_err = 0.0;
   double sum_err_sqrd = 0.0;
   double delta_phase = (2.0 * PI) / num_points;
   for (ii = 0; ii < num_points; ii++) {
      phase = delta_phase * ii;
      real = cos(phase);
      imag = sin(phase);
      err = sqrt(real * real + imag * imag)
            - alpha_beta_mag(alpha, beta, real, imag);
      sum_err += err;
      sum_err_sqrd += err * err;
      err = fabs(err);
      if (err > peak_err) {
         peak_err = err;
      }
   }
   avg_err = sum_err / num_points;
   rms_err = sqrt(sum_err_sqrd / num_points);
   printf("%-16s %14.12lf %14.12lf  %9.6lf %4.1lf %4.1lf\n",
          name, alpha, beta, avg_err, decibels(rms_err),
          decibels(peak_err));
}
/*********************************************************************/
void main(void)
{
   #define NUM_CHECK_POINTS 100000
   typedef struct tagALPHA_BETA {
      char *name;
      double alpha;
      double beta;
   } ALPHA_BETA;
   #define NUM_ALPHA_BETA 16
   const ALPHA_BETA coeff[NUM_ALPHA_BETA] = {
      { "Min RMS Err",      0.947543636291, 0.3924854250920 },
      { "Min Peak Err",     0.960433870103, 0.3978247347593 },
      { "Min RMS w/ Avg=0", 0.948059448969, 0.3926990816987 }, 
      { "1, Min RMS Err",              1.0,     0.323260990 },
      { "1, Min Peak Err",             1.0,     0.335982538 },
      { "1, 1/2",                      1.0,      1.0 / 2.0  },
      { "1, 1/4",                      1.0,      1.0 / 4.0  },
      { "Frerking",                    1.0,            0.4  },
      { "1, 11/32",                    1.0,     11.0 / 32.0 },
      { "1, 3/8",                      1.0,      3.0 / 8.0  },
      { "15/16, 15/32",        15.0 / 16.0,     15.0 / 32.0 },
      { "15/16, 1/2",          15.0 / 16.0,      1.0 / 2.0  },
      { "31/32, 11/32",        31.0 / 32.0,     11.0 / 32.0 },
      { "31/32, 3/8",          31.0 / 32.0,      3.0 / 8.0  },
      { "61/64, 3/8",          61.0 / 64.0,      3.0 / 8.0  },
      { "61/64, 13/32",        61.0 / 64.0,     13.0 / 32.0 }
   };
   int ii;
   printf("\n             Alpha * Max + Beta * Min Magnitude Estimator\n\n");
   printf("Name                  Alpha           Beta       Avg Err  RMS   Peak\n");
   printf("                                                 (linear) (dB)  (dB)\n");
  
printf("---------------------------------------------------------------------\n");
   for (ii = 0; ii < NUM_ALPHA_BETA; ii++) {
      test_alpha_beta(coeff[ii].name, coeff[ii].alpha, coeff[ii].beta, 1024);
   }
}

Best Answer

Try this one - it is not debugged and written in the forum's interface:

function main
   NUM_CHECK_POINTS = 100000;
   coeff = { ...
       'Min RMS Err',      0.947543636291, 0.3924854250920; ...
       'Min Peak Err',     0.960433870103, 0.3978247347593; ...
       'Min RMS w/ Avg=0', 0.948059448969, 0.3926990816987; ...
       '1, Min RMS Err',              1.0,     0.323260990; ...
       '1, Min Peak Err',             1.0,     0.335982538; ...
       '1, 1/2',                      1.0,      1.0 / 2.0 ; ...
       '1, 1/4',                      1.0,      1.0 / 4.0 ; ...
       'Frerking',                    1.0,            0.4 ; ...
       '1, 11/32',                    1.0,     11.0 / 32.0; ...
       '1, 3/8',                      1.0,      3.0 / 8.0 ; ...
       '15/16, 15/32',        15.0 / 16.0,     15.0 / 32.0; ...
       '15/16, 1/2',          15.0 / 16.0,      1.0 / 2.0 ; ...
       '31/32, 11/32',        31.0 / 32.0,     11.0 / 32.0; ...
       '31/32, 3/8',          31.0 / 32.0,      3.0 / 8.0 ; ...
       '61/64, 3/8',          61.0 / 64.0,      3.0 / 8.0 ; ...
       '61/64, 13/32',        61.0 / 64.0,     13.0 / 32.0};
   
   fprintf('\n             Alpha * Max + Beta * Min Magnitude Estimator\n\n');
   fprintf('Name                  Alpha           Beta       Avg Err  RMS   Peak\n');
   fprintf('                                                 (linear) (dB)  (dB)\n');
  
printf('---------------------------------------------------------------------\n');
   for ii = 1:size(coeff, 1)
      test_alpha_beta(coeff{ii, 1}, coeff{ii, 2}, coeff{ii, 3}, 1024);
   end
end
function R = alpha_beta_mag(alpha, beta, inphase, quadrature)
% magnitude ~= alpha * max(|I|, |Q|) + beta * min(|I|, |Q|) */
   abs_inphase    = abs(inphase);
   abs_quadrature = abs(quadrature);
   if abs_inphase > abs_quadrature
      R = alpha * abs_inphase + beta * abs_quadrature;
   else
      R = alpha * abs_quadrature + beta * abs_inphase;
   end
end
function R = decibels(linear)
   R = 20.0 * log10(max(1e-20, linear));
end
function R = test_alpha_beta(name, alpha, beta, num_points)
   peak_err = 0.0;
   sum_err = 0.0;
   sum_err_sqrd = 0.0;
   delta_phase = (2.0 * Ppi) / num_points;
   for ii = 0:num_points - 1
      phase = delta_phase * ii;
      real = cos(phase);
      imag = sin(phase);
      err = sqrt(real * real + imag * imag) ...
            - alpha_beta_mag(alpha, beta, real, imag);
      sum_err = sum_err + err;
      sum_err_sqrd = sum_err_sqrd + err * err;
      err = abs(err);
      if (err > peak_err)
         peak_err = err;
      end
   end
   
   avg_err = sum_err / num_points;
   rms_err = sqrt(sum_err_sqrd / num_points);
   fprintf('%-16s %14.12f %14.12f  %9.6f %4.1f %4.1f\n',
          name, alpha, beta, avg_err, decibels(rms_err),
          decibels(peak_err));
end

It took 2 minutes to convert is from C to Matlab.

Related Solutions

MATLAB: Can someone rearrange the code to run

%subdivisions space /time
ht=0.01; Tmax=1.2; nx=33; hx=1/(nx-1); x=[0:hx:1]';
%matrices
K=stiff2(1/pi^2,hx,nx); M=mass(1/ht,hx,nx); 
A=M+K;
%boundary conditions by deleting rows
A(nx,:)=[]; A(:,nx)=[]; A(1,:)=[]; A(:,1)=[];
%creation
R=chol(A);
%initial condition
u=cos(pi*(x-1/2));
%time step loop
k=0; 
while (k*ht<Tmax)
k=k+1;
b=M*u;
b(nx)=[];b(1)=[];
%solve
u=zeros(nx,1);
u(2:nx-1)=R\(R'\b);
end
% Set the matrices to zero
k=zeros(6,6,2); K=zeros(6,6,2); Gamma=zeros(6,6,2);
% Enter parameter values, in units: mm^2, mm^4, and MPa(10^6 N/m)
A=6000; II=200*10^6; EE=200000;
% Convert units into meters and kN
A=A/10^6; II = II/10^12; EE =EE*1000;
% Element 1
i=[0,0]; j=[7.416,3];
[k(:,:,1),K(:,:,1),Gamma(:,:,1)]=stiff(EE,II,A,i,j);
% Element 2
i=j; j=[15.416,3];
[k(:,:,2),K(:,:,2),Gamma(:,:,2)]=stiff(EE,II,A,i,j);
% Define the elementary stiffness matrix
ID=[-4 1 -7;-5 2 -8;-6 3 -9];
% Define the connection matrix
LM=[-4 -5 -6 1 2 3; 1 2 3 -7 -8 -9];
% Assemble the augmented stiffness matrix
Kaug = zeros(9); 
for elem=1:2
for r=1:6
lr=abs(LM(elem,r));
for c=1:6
lc=abs(LM(elem,c));
Kaug(lr,lc)=Kaug(lr,lc)+K(r,c,elem);
end
end
end
% Extract the stiffness matrix
Ktt=Kaug(1:3,1:3);
% Determine the reactions at the nodes in local coordinates
fea(1:6,1)=0; 
fea(1:6,2)=[0,8*4/2,4*8^2/12,0,8*4/2,-4*8^2/12]';
% Determine the reactions in the global coordinate system
FEA(1:6,1)=Gamma(:,:,1)'*fea(1:6,1);
FEA(1:6,2)=Gamma(:,:,2)'*fea(1:6,2);
% FEA_Rest for constrained nodes
FEA_Rest=[0,0,0,FEA(4:6,2)'];
% Assemble the right-hand side for non-constrained nodes
P(1)=50*3/8; P(2)=-50*7.416/8-FEA(2,2); P(3)=-FEA(3,2);
% Solve to find the displacements in meters and in radians
Displacements=inv(Ktt)*P'
% Extract Kut
Kut = Kaug(4:9,1:3);
% Compute the reactions and introduce boundary conditions
Reactions=Kut*Displacements+FEA_Rest'
% Solve to find the internal forces, excluding fixed points
dis_global(:,:,1)=[0,0,0,Displacements(1:3)'];
dis_global(:,:,2)=[Displacements(1:3)',0,0,0]; 
for elem=1:2
dis_local= Gamma(:,:,elem)*dis_global(:,:,elem)';
int_forces= k(:,:,elem)*dis_local+fea(1:6,elem)
end
%The above script calls the stiff function, which can be implemented as follows:
function [k,K,Gamma] = stiff( EE,II,A,i,j )
% Find the length
L=sqrt((j(2)-i(2))^2+(j(1)-i(1))^2);
% Compute the angle theta
if(j(1)-i(1))~=2
alpha=atan((j(2)-i(2))/(j(1)-i(1)))
else
alpha=-pi/2;
end
% Form the rotation matrix Gamma
Gamma =[cos(alpha) sin(alpha) 0 0 0 0;
-sin(alpha) cos(alpha) 0 0 0 0;
0 0 1 0 0 0;
0 0 0 cos(alpha) sin(alpha) 0;
0 0 0 -sin(alpha) cos(alpha) 0;
0 0 0 0 0 1];
% Form the elementary stiffness matrix in local coordinates
EI=EE*II; EA=EE*A; 
k=[EA/L, 0, 0, -EA/L, 0, 0; 
    0, 12*EI/L^3, 6*EI/L^2, 0, -12*EI/L^3,6*EI/L^2;
0, 6*EI/L^2, 4*EI/L, 0 -6*EI/L^2, 2*EI/L;
-EA/L, 0 ,0 , EA/L, 0, 0;
0, -12*EI/L^3, -6*EI/L^2, 0, 12*EI/L^3, -6*EI/L^2;
0, 6*EI/L^2, 2*EI/L, 0, -6*EI/L^2, 4*EI/L];
% Elementary matrix in global coordinates
K=Gamma'*k*Gamma;
end
function R=stiff2(nu,h,n)
%
[m1,m2]=size(nu);
if (length(nu)==1)
ee=nu*ones(n-1,1)/h; e1=[ee;0]; e2=[0;ee]; e=e1+e2;
else
ee=.5*(nu(1:n-1)+nu(2:n))/h; e1=[ee;0]; e2=[0;ee]; e=e1+e2;
end
R=spdiags([-e1 e -e2],-1:1,n,n);
return
end
function M=mass(alpha,h,n)
[m1,m2]=size(alpha);
if(length(alpha)==1)
ee=alpha*h*ones(n-1,1)/6;
e1=[ee;0]; e2=[0;ee]; e=e1+e2;
else
ee=h*(alpha(1:n-1)+alpha(2:n))/12;
e1=[ee;0]; e2=[0;ee]; e=e1+e2;
end
M=spdiags([e1 2*e e2],-1:1,n,n);
return
end

MATLAB: There is the error that the number of input parameters is insufficient

integral2(Tuu,0.01,pi/2,0,pi/4)

is equivalent for this purpose to

integral2(Tuu(),0.01,pi/2,0,pi/4)

which is to say that Tuu is called with no parameters, and is expected to output a function handle that will then be integrated with integral2.

You need

integral2(@Tuu,0.01,pi/2,0,pi/4)

Best Answer

Related Solutions

MATLAB: Can someone rearrange the code to run

MATLAB: There is the error that the number of input parameters is insufficient

Related Question