MATLAB Answers

Subtracting complex number from a discrete signal

10 views (last 30 days)
Michelle Weinmann on 11 May 2021
Edited: Michelle Weinmann on 20 May 2021 at 15:16
Hello Matlab community,
I'm trying to remove a background signal from a set of discrete signals. I have determined that the background signal is represented by -0.000053 + 0.0030995i, or 0.0031*e^(i*1.588). In my code, the "j" loop represents each of the 11 positions where signals were recorded. I'm having trouble understanding how to actually perform the subtraction of the background signal from the 11 discrete signals; it seems that I need to convert the discrete digital signals into the complex plane but I don't know how to do that. Any assistance would be appreciated; I've attached the code.
Thank you,
Michelle
% DataF.mat has 250000 data points per position. Each point of data
% corresponds to 1mm of motion beginning at -5mm from center to 5mm from
% center ranging 10mm total, resulting in 11 data points.
%
clc;
clear;
close all;
%% Constants
EPSfreq=2.36e3; %EPS carrier
rate=50000; %LabView data file Sample Clock rate = rate
dt=1/rate;
points=250000;
%% Main Code
load('DataF.mat');
% Initialized Matrices
% [c,t,tm,out]=deal(zeros(points,numPoints));
% nSamp=Temp_points
for i=1:points
t(i)=i*dt;
tm(i)=t(i)*1e6; %Time in microseconds
end
%Position loop: 11 positions, 100 points per measured location
for j=1:11
%Output from filter (pickup signal)
out(:,j) = DataF((j-1)*points+1:(j*(points)), 1);
%Reference from excitation coils
c(:,j) = DataF((j-1)*points+1:(j*(points)), 2);
%Function detrend removes data bias by computing mean values and subtrancting them from signals
%Use if DAQ system puts offset on signals
out_d(:,j) = detrend(out(:,j));
c_d(:,j) = detrend(c(:,j));
%Plot of detrended signals: plots for each position
% figure
% plot(t,c_d(:,j),t,out_d(:,j));
% xlabel('Time (s)')
% ylabel('Amplitude (V)')
%Cross-correlation and auto-correlation
%Auto correlation reference
Rx(j) = (1/points)*sum(c_d(:,j).^2);
%Auto correlation measured
Ry(j) = (1/points)*sum(out_d(:,j).^2);
%Amplitude of reference
Ax(j) = sqrt(2*Rx(j));
%Amplitude of measured
Ay(j) = sqrt(2*Ry(j));
%Cross-correlation
R(j) = (1/points)*sum(c_d(:,j).*out_d(:,j));
%Relative phase in degrees
relphase(j)=acos(R(j)/sqrt(Rx(j)*Ry(j)))*(180/pi);
%Amplitude ratio
relamp(j) = Ay(j)/Ax(j);
end
%Sign adjust loop
for k=1:5 %First 5 positions
signadj(k)=sign(relphase(k)-180); %Adjusts amplitude graph slope
relamp(k)=relamp(k)*signadj(k); %Output results in linear graph slope
end
figure
plot([-5:1:5],relamp,'-*');grid
xlabel('Displacement')
ylabel('Derived relative amplitude')
% axis([-5 5 -0.05 0.05])
figure
plot([-5:1:5],relphase,'-*');grid
xlabel('Displacement')
ylabel('Phase difference')
figure
subplot(2,2,1);
plot(Rx)
title('Auto Correlation Sum: Excitation')
legend('X')
grid on
subplot(2,2,2);
hold on
plot(R)
plot(Ry)
title('Auto and Cross Correlation Sum: Cross correlation and Pickup')
legend('R','Y')
figure
subplot(2,2,1);
plot(Ax)
title('Excitation Amplitude')
legend('Ax')
grid on
subplot(2,2,2);
plot(Ay)
title('Pickup Amplitude')
legend('Ay')
grid on
9 CommentsShowHide 8 older comments
Mathieu NOE on 13 May 2021
if you have electric perturbations from mains , this happen to known and fix frequencies (50 or 60 Hz and maybe a bit of harmonics)
why not simply use a notch filter and do the filtering directly on the time data ?
example attached

Sign in to comment.

Answers (2)

Walter Roberson on 12 May 2021
A = 123
A = 123
B = -0.000053 + 0.0030995i
B = -0.0001 + 0.0031i
C = A - B
C = 1.2300e+02 - 3.0995e-03i
In other words, you do not need to do anything special to convert the signal into the complex plane: MATLAB will automatically extend real values to complex as needed.
1 CommentShowHide None
Walter Roberson on 12 May 2021
You seem to be under the impression that you must do something specific to convert your discrete signal to complex representation. If so then
complex_signal = complex(real_signal, 0)
If you need complex form A*exp(j*theta) then this is the same as theta = 0, as clearly exp(j*0) --> exp(0) -> 1 and A*1 is A
However, as I pointed out, in MATLAB you can simply subtract the complex value and MATLAB will automatically extend the real value to complex .

Sign in to comment.

Michelle Weinmann on 20 May 2021 at 15:07
Edited: Michelle Weinmann on 20 May 2021 at 15:16
Hello all,
It turns out I had the wrong idea about what I needed to do. I needed to subtract a complex number representing a signal from another complex number representing the output signal using the output values "relamp" as amplitude and "relphase" as phase. I have adjusted my code as required and attached it for reference.
However, it's not performing quite how it should. According to my advisor, the "relamp" values should be quite close to linear and I have redundant steps in the subtraction step which changes some of the values slightly; the linear relationship is slightly off. I can't see where the redundant steps are and how they would affect the relamp values. I've also attached a shortened data file (DataResized.mat) and the plot resulting (RelativeAmplitudePlot.fig); the first half of the points look ok but the second half do not. I've also attached a photo of what the plot should look like (NewPlot.jpg). Any help would be appreciated.
Thank you,
Michelle
% DataF.mat has 250000 data points per position. Each point of data
% corresponds to 1mm of motion beginning at -5mm from center to 5mm from
% center ranging 10mm total, resulting in 11 data points.
%
clc;
clear;
close all;
%% Constants
EPSfreq=2.36e3; %EPS carrier
rate=50000; %LabView data file Sample Clock rate = rate
dt=1/rate;
points=2500;
%% Main Code
load('DataResized.mat');
% Initialized Matrices
% [c,t,tm,out]=deal(zeros(points,numPoints));
% nSamp=Temp_points
for i=1:points
t(i)=i*dt;
tm(i)=t(i)*1e6; %Time in microseconds
end
%Position loop: 11 positions, 100 points per measured location
for j=1:11
%Output from filter (pickup signal)
out(:,j) = DataRsz((j-1)*points+1:(j*(points)), 1);
%Reference from excitation coils
c(:,j) = DataRsz((j-1)*points+1:(j*(points)), 2);
%Function detrend removes data bias by computing mean values and subtrancting them from signals
%Use if DAQ system puts offset on signals
out_d(:,j) = detrend(out(:,j));
c_d(:,j) = detrend(c(:,j));
%Plot of detrended signals: plots for each position
% figure
% plot(t,c_d(:,j),t,out_d(:,j));
% xlabel('Time (s)')
% ylabel('Amplitude (V)')
%Cross-correlation and auto-correlation
%Auto correlation reference
Rx(j) = (1/points)*sum(c_d(:,j).^2);
%Auto correlation measured
Ry(j) = (1/points)*sum(out_d(:,j).^2);
%Amplitude of reference
Ax(j) = sqrt(2*Rx(j));
%Amplitude of measured
Ay(j) = sqrt(2*Ry(j));
%Cross-correlation
R(j) = (1/points)*sum(c_d(:,j).*out_d(:,j));
%Correlation ratio for relphase
Rrat(j) = R(j)/sqrt(Rx(j)*Ry(j));
%Relative phase gives angles between 0 and 180 deg
relphase(j)=acos(Rrat(j));
%Amplitude ratio
relamp(j) = Ay(j)/Ax(j);
%Subtract background signal in complex plane
V1(j)= relamp(j)*exp(1i*(relphase(j)));
%Background signal
V2 = -0.00005332875 + 0.0030995i;
V3(j) = V1(j)-V2;
%Update relamp and relphase values once background signal subtracted
relamp(j) = abs(V3(j));
relphase(j) = angle(V3(j))*180/pi;
%Update pickup amplitude once background signal subtracted
Ay(j) = relamp(j)*Ax(j);
end
%Sign adjust loop: since acos only outputs angles from 0 to 180
for k=1:11
if Rrat(k) > 0
relphase(k) = relphase(k) - 2*pi;
end
if Rrat(k) < 0
relamp(k) = -relamp(k);
end
end
figure
plot([-5:1:5],relamp,'o-');grid
legend('Adjusted pickup relative amplitude','V3')
xlabel('Displacement')
ylabel('Derived relative amplitude')
0 CommentsShowHide -1 older comments

Sign in to comment.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!