Measuring average velocity with light gate signal

1 次查看（过去 30 天）
Elon Hendriksen
Elon Hendriksen 2022-10-10
回答： TsungYeh Chou 2022-11-13
Hi all!
I have a simple issue which i cannot solve for some reason and I hope someone can help me here.
I have a signal from a light gate which returns a noisy signal between 4 and 5 volts when the light gate is uninterupted. When something blocks the gate the signal goes to 0V.
I would like to find the width of the interuptions to calculate the average velocity of the object interupting the light gate signal.
I have tried to use the [pks,locs,w,p] = findpeaks(data); function but it does not seem to work
clear
close all
clc
%% import data
%fileName = "Flatplate.tdms";
fileName = "Wig15deg.tdms";
%fileName = "Wig30deg.tdms";
info = tdmsinfo(fileName);
info.ChannelList;
tdmsreadprop(fileName);
group = "Untitled";
channel = "Pres_1";
tdmsreadprop(fileName, ChannelGroupName=group, ChannelName=channel);
time = 10 % [s]
freq = 100000 % [Hz]
volt = 5
maxpres = 68.95 % [Bar]
maxacc = 500 % [g]
time_output = [0.0, 1] %, 10]
dist_light = 0.4
val_light = 2.0
g = -9.81
%% Data handling
timeStep = tdmsreadprop(fileName, ChannelGroupName=group, ChannelName=channel, PropertyNames="wf_increment");
startTime = tdmsreadprop(fileName, ChannelGroupName=group, ChannelName=channel, PropertyNames="wf_start_time");
data = tdmsread(fileName, StartTime=startTime.wf_start_time, TimeStep=seconds(timeStep.wf_increment));
ttData = data{1};
Double = table2array(ttData);
% cut in to separate collums
Pres1=Double(:,2)/ volt * maxpres * 1e5;
Pres2=Double(:,3)/ volt * maxpres * 1e5;
Pres3=Double(:,4)/ volt * maxpres * 1e5;
Pres4=Double(:,5)/ volt * maxpres * 1e5;
Acc1=Double(:,6)/ volt * maxacc;
Acc2=Double(:,7)/ volt * maxacc;
%find max and cut columm into section before and after impact
maximum = max(max(Pres1));
[P1x]=find(Pres1==maximum);
maximum = max(max(Pres2));
[P2x]=find(Pres2==maximum);
maximum = max(max(Pres3));
[P3x]=find(Pres3==maximum);
maximum = max(max(Pres4));
[P4x]=find(Pres4==maximum);
peak_min = 100;
peak_plus = 300;
CleanPres1=Pres1(P1x-peak_min:P1x+peak_plus,1);
CleanPres2=Pres2(P2x-peak_min:P2x+peak_plus,1);
CleanPres3=Pres3(P2x-peak_min:P2x+peak_plus,1);
CleanPres4=Pres4(P1x-peak_min:P1x+peak_plus,1);
CleanAcc1=Acc1(P1x-peak_min:P1x+peak_plus,1);
%all data back into collumns
P_all=[CleanPres1 CleanPres2 CleanPres3 CleanPres4];
P_all2=[CleanPres1 CleanPres4];
%% Symmetry check
symmetry = rmoutliers((CleanPres3)./(CleanPres2));
avg_sym_inside = mean(symmetry, 'all')
symmetry2 = rmoutliers((CleanPres1)./(CleanPres4));
avg_sym_outside = mean(symmetry2, 'all')
%% Velocity measurement
%signal = [0 0 6 6 6 0 0 6 6 6 0 0 6 6 6 0 0]
signal=Double(:,8);
halfMaxValue = 0.5 * (min(signal) + max(signal)); % should be around 2.5
lowSignals = signal < halfMaxValue;
diffSignal = diff(lowSignals);
[pks,locs,w,p] = findpeaks(diffSignal);
% diffSignalIndexes = find(diffSignal == 0);
% [~,peaklocs] = findpeaks(diffSignalIndexes);
% tstep = 1/freq;
% period = (diff(peaklocs))*tstep
% diffSignalIndexes = find(diffSignal == 0);
% period = diffSignalIndexes(2) - diffSignalIndexes(1) % Distance between first two pulses.
%% Show data
figure(1)
stackedplot(P_all)
figure(2)
stackedplot(ttData)
figure(3)
plot(symmetry)
figure(4)
plot(symmetry2)
figure(5)
plot(P_all2)
figure(6)
plot(signal)
Below a screenshot of my data:

请先登录，再回答此问题。

回答（2 个）

Image Analyst
Image Analyst 2022-10-10
Ran in:
Unfortunately you forgot to attach your data so I'll have to guess. If you have the Image Processing Toolbox, try this:
signal = rand(1, 50); % Sample data
% Find elements where the signal is less than 0.5.
mask = signal < 0.5;
% Measure the widths of those signals
props = regionprops(mask, 'Area')
props = 11×1 struct array with fields:
Area
allInterruptionWidths = [props.Area]
allInterruptionWidths = 1×11
1 2 1 1 1 1 1 1 3 6 1
meanWidth = mean(allInterruptionWidths)
meanWidth = 1.7273
If you have any more questions, then attach your data and code to read it in with the paperclip icon after you read this:
TUTORIAL: How to ask a question (on Answers) and get a fast answer
  1 个评论
Elon Hendriksen
Elon Hendriksen 2022-10-11
编辑：Elon Hendriksen 2022-10-11
Ah I'm sorry, I have added my code now! Unfortunatley i cant send the TDMS file since it exceeds 5MB, even after zipping. I have added an excel file with half of my lightgate data.
I checked your code proposal and it seems to work! thanks!

请先登录，再进行评论。

TsungYeh Chou
TsungYeh Chou 2022-11-13
Hello Elon,
May I ask which light gate (brand and model) do you use? And how do you connect the light gate to matlab? I want to use lightgate to monitor a event and have matlab to perfrom a script. Thank you!

类别

在 Help Center 和 File Exchange 中查找有关 Annotations 的更多信息

标签

产品


版本

R2022a

Community Treasure Hunt

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

Start Hunting!

Translated by