My code keeps telling me with an error "Index Exceeds Array Bounds" HELP PLEASE

3 visualizaciones (últimos 30 días)
Mohammad Junayed
Mohammad Junayed el 25 de Oct. de 2018
Respondida: Mohammad Junayed el 25 de Oct. de 2018
%% Initialization clc clear close all
%% Data Def AppOpenCnt = 0; PreEventCnt = 0; PreE = cell(0);%Pre event EventTkt = 0; LC_flag = 0; d = cell(0); % Data record e = cell(0);% Event aa = cell(0);% Applicance database
nFreqDAQ = 9; gMinThr = 30;%Minimum power trigger threshold gMinRate = 0.2;%Power change trigger ratio gStdP = 5;%%Step size of normal power flactuation gStdI = 1;%%Step size of normal amplitude fluctuation %Normalization coefficient wP = 0.9; wQ = 0; wIh = 0.1; nH = floor(nFreqDAQ/2); ASDPath = 'user.xlsx'; %DataPath = '60V20AH.txt'; %%DataPath = '36V12AH.txt'; %%DataPath = '48V12AH.txt'; DataPath = 'data\20181015.xls'; % DataPath = 'test0930.xls'; PrePara.PreThreshold = 5;%Pre-determination window PrePara.MinThr = 10;%Pre-judgment trigger threshold, can be judged PrePara.prejudgeFlag = 0; PrePara.PreSteaTkt = 0; PrePara.PreEventCnt = 0; PrePara.PrejudgeTkt = 0; PrePara.PreE = cell(0); Threshold = 5;%Decision window unit seconds
dPRateThr = 0.7;%%Microwave oven standby avg = 23w, the maximum deviation is 0.267 GAvgdP = 96; GavgdQ = 40; GavgIh = 10;
bText = 1; %% Read ASD File IsInRecord = false; iState = -1; iLogic = -1; f = fopen(ASDPath, 'r'); idx = 1; while ~feof(f) s = fgetl(f); if isempty(s) continue; end if length(s) < 5 continue; end if s(1) == '%' continue; end if ~IsInRecord % Read Record Header ss = strtrim(split(s, ',')); tt.Type = ss{1}; tt.Brand = str2double(ss{2}); tt.Model = str2double(ss{3}); tt.Note = ss{4}; tt.RateW = str2double(ss{5}); tt.IsPlural = (ss{6} == 1); tt.nFreq = str2double(ss{7}); tt.nState = str2double(ss{8}); tt.IsSwtichable = (ss{9} == 1); tt.nLogic = str2double(ss{10}); tt.State = cell(tt.nState,1); tt.Logic = cell(tt.nLogic,1); IsInRecord = true; iState = 0; else % Read Record Content if iState >= 0 && iState < tt.nState ss = strtrim(split(s, ',')); tt.StateName = ss{1}; tt.V = str2double(ss{2}); tt.Im = zeros(ceil(tt.nFreq/2),1); tt.Ia = zeros(ceil(tt.nFreq/2),1); tt.Ih = zeros(floor(tt.nFreq/2),1); tt.I = zeros(floor(tt.nFreq/2),1); % tt.IHD = zeros(floor(t.nFreq/2),1); for i = 1:ceil(tt.nFreq/2) tt.Im(i) = str2double(ss{i*2+1}); tt.Ia(i) = str2double(ss{i*2+2}); if i > 1 tt.Ih(i-1) = tt.Im(i); end end tt.P = tt.V*tt.Im(1)*cos(tt.Ia(1)); tt.Q = tt.V*tt.Im(1)*sin(tt.Ia(1)); iState = iState + 1; tt.State{iState} = tt; if iState == tt.nState iLogic = 0; end elseif iLogic >= 0 && iLogic < tt.nLogic error('Special logic not implemented'); iLogic = iLogic + 1; else error('Status code error'); end % IsFinished if iState == tt.nState && iLogic == tt.nLogic && IsInRecord aa{end+1,1} = tt; %#ok<SAGROW> iState = -1; iLogic = -1; IsInRecord = false; clear tt end end end fclose(f); fprintf('Read ASD completed, %d records imported.\n', size(aa,1)); r = cell(size(aa,1),1);
%% Event visualition file = 'red.jpg'; pic = imread(file); pic_gray = rgb2gray(pic); pic_green = imread('green.jpg'); figAll = figure('numbertitle','off','name','Event Detection');figAll.Position=[700 100 1200 800]; subplot(4,3,1);figLED = 1 ; imshow(pic_gray);title('LED','FontSize',50);hold on; subplot(4,3,2);figCFL = 2 ; imshow(pic_gray);title('CFL','FontSize',50);hold on; subplot(4,3,3);figWave = 3 ; imshow(pic_gray);title('Microwave','FontSize',50);hold on; subplot(4,3,4);figAuto = 4 ; imshow(pic_gray);title('Iron','FontSize',50);hold on; subplot(4,3,5);figKettle = 5 ; imshow(pic_gray);title('Kettle','FontSize',50);hold on; subplot(4,3,6);figHairdryer = 6 ; imshow(pic_gray);title('Hairdryer','FontSize',50);hold on; subplot(4,3,7);figCooker = 7 ; imshow(pic_gray);title('Cooker','FontSize',50);hold on; subplot(4,3,8);figHeater = 8 ; imshow(pic_gray);title('Heater','FontSize',50);hold on; subplot(4,3,9);figBath = 9 ; imshow(pic_gray);title('Bath','FontSize',50);hold on; subplot(4,3,11);figLeakCurrent = 11 ; imshow(pic_gray);title('LeakCurrent','FontSize',50);hold on; drawnow;
%% Data Loading %%Message format Time Frequency Fundamental voltage amplitude Leakage current Fundamental /3/5/7/9th harmonic current amplitude Fundamental /3/5/7/9th harmonic phase angle %% 2018-10-02 14:50:48 50.0355 233.2042 34.8354 0.2678 0.2589 0.2185 0.1704 0.1244 25.1995 -83.8638 108.7787 71.9556 -24.6973 ,- P = 56.509, Q = -26.590 [a,b,c]=fileparts(DataPath); qq = lower(c); if(c == ".txt") bText = 1; else bText = 0; end
if(bText) s = importdata(DataPath); ss = s.textdata(:,2:20); else
[TXT, NUM] = xlsread(DataPath); ss = NUM; end %%% para Initialize [PrePara,EventTkt] = ParaInit(PrePara,EventTkt);
for i = 1:length(ss) if(bText) tt.Time = datetime(str2double(ss{i,1}),str2double(ss{i,2}),str2double(ss{i,3}),str2double(ss{i,4}),str2double(ss{i,5}),str2double(ss{i,6})); tt.Freq = str2double(ss{i,7}); tt.V = str2double(ss{i,8}); tt.Ir = str2double(ss{i,9});%Residual current RMS value tt.Im = str2double(ss{i,10});%Fundamental current amplitude tt.Ia = str2double(ss{i,15});%Fundamental current phase angle tt.Ih = zeros(nH,1);%%Harmonic amplitude tt.Iha = zeros(nH,1);%%Harmonic phase angle tt.I = zeros(nH,1);%%Vector difference from the target value harmonic tt.P = tt.V*tt.Im*cos((tt.Ia/180)*pi); %active tt.Q = tt.V*tt.Im*sin((tt.Ia/180)*pi); %reactive fprintf(' - P = %.3f, Q = %.3f, Current Data is %d \n', tt.P, tt.Q, i); for iH = 1:nH tt.Ih(iH) = abs(str2double(ss{i,10+iH})); tt.Iha(iH) = abs(str2double(ss{i,15+iH})); end else tt.Freq = 50.00; tt.V = 220.00; tt.Im = str2double(NUM{i, 4}); tt.Ia = str2double(NUM{i, 9}); %Fundamental current phase angle tt.Ir = 0; %Residual current RMS value tt.Ih = zeros(nH,1); %Harmonic power tt.Iha = zeros(nH,1); %%Harmonic phase angle tt.I = zeros(nH,1); %%Vector difference from the target value harmonic tt.P = str2double(NUM{i, 2}); %Active tt.Q = str2double(NUM{i, 3}); %Reactive tt.Time = NUM(i,14); fprintf(' - P = %.3f, Q = %.3f, Ia = %.3f, Current Data is %d \n', tt.P, tt.Q, tt.Ia, i);
for iH = 1:nH
tt.Ih(iH) = abs(str2double(NUM{i,4+iH}));
end
if (abs(tt.Freq - 50)>=1) || (abs(tt.V-220)>=22)
fprintf(' - System not ready, wait for the next message\n');
clear tt
% continue; % Noise, skip the following procedure
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
d{end+1,1} = tt;
clear tt
% Edge ==> Event
IsEventNow = false;
%%Leakage current detection
LC_flag = IrDect(figAll,figLeakCurrent,pic,pic_gray,d);
if size(d,1) > 2
dP = d{end}.P - d{end-1}.P;%Active difference
dP2 = d{end-1}.P - d{end-2}.P;
dQ = d{end}.Q - d{end-1}.Q;
dIh = abs(d{end}.Ih - d{end-1}.Ih);
dIr = d{end}.Ir - d{end-1}.Ir;
%%%%%%%%%%%%%Debugging
if i ==196
disp('Open time!');
end
if i ==508
disp('Open time!');
end
if i ==102
disp('Open time!');
end
%%%%%%%%%%%%%Prejudging begins
if PrePara.prejudgeFlag == 0
dPRate = abs((abs(dP)-abs(dP2))/dP) ;
if (dPRate > dPRateThr) && (abs(max(abs(dP),abs(dP2))) > PrePara.MinThr)
fprintf(' - Edge detected: dP = %.3f, dQ = %.3f, ||dIh|| = %.3f\n', dP, dQ, norm(dIh));
PrePara.prejudgeFlag = 1;
PrePara.PrejudgeTkt = size(d,1) -1;
else
% LC_flag = IrDect(figAll,figLeakCurrent,pic,pic_gray,d);
end
elseif PrePara.prejudgeFlag == 1%%Confirmation of judgment
PrePara.PreEventCnt = PrePara.PreEventCnt + 1;
fprintf('Event continue No.: %.3f',PrePara.PreEventCnt);
% t.Start = EventTkt;
% t.End = size(d,1); % two-end inclusive
tt.dP = d{end}.P - d{PrePara.PrejudgeTkt}.P;
tt.dQ = d{end}.Q - d{PrePara.PrejudgeTkt}.Q;
tt.dIr = d{end}.Ir - d{PrePara.PrejudgeTkt}.Ir;
tt.dIh = abs(d{end}.Ih - d{PrePara.PrejudgeTkt}.Ih);
for h = 1:nH
[tt] = CalHarmVec(tt,d,h,PrePara.PrejudgeTkt);
end
PrePara.PreE{end+1,1} = tt;
PrePara.PreErrCnt = 0;
PreRate = 0;
%%calculate the average value of the latest 6 seconds
if(PrePara.PreEventCnt > PrePara.PreThreshold)
[AvgdP,AvgdQ,AvgdIh] = AvgCalc(PrePara.PreEventCnt,PrePara.PreE);
%%Judge steady state coming or not
for j = 1:PrePara.PreEventCnt
PrePara.PreSteaTkt = PrePara.PrejudgeTkt + j + 1;
PredP = abs(AvgdP);
PreRate = abs(PrePara.PreE{end-PrePara.PreEventCnt + PrePara.PreErrCnt +1,1}.dP) /PredP;
%%Entering steady state
if (abs(1 - PreRate) < gMinRate) && (PredP > gMinThr)
%%If the data sample is not enough, jump out
if j > ceil(PrePara.PreThreshold/2)
[PrePara,EventTkt] = ParaInit(PrePara,EventTkt);
break;
end
if EventTkt == 0
EventTkt = PrePara.PrejudgeTkt;
IsEventNow = true;
% fprintf(' - Event created at %s\n', d{end,1}.Time);
e = PrePara.PreE;
[PrePara,EventTkt] = ParaInit(PrePara,EventTkt);
break;
else
fprintf(' - Event in progress at %s\n', d{end,1}.Time);
end
else
[AvgdP,AvgdQ,AvgdIh] = AvgCalc((PrePara.PreEventCnt - j),PrePara.PreE);
PrePara.PreErrCnt = PrePara.PreErrCnt + 1;
if j == PrePara.PreEventCnt
[PrePara,EventTkt] = ParaInit(PrePara,EventTkt);
end
end
end
end
%%Leakage current detection
% LC_flag = IrDect(figAll,figLeakCurrent,pic,pic_gray,d);
if i ==34
disp('Open time!');
end
if i ==508
disp('Close time!');
end
% Event Process
if IsEventNow
selObj = -Inf;
selA = 0; % Best fit appliance
selS = 0; % Best fit state
selP = zeros(1,3);
Sign = 2*(e{end}.dP > 0)-1;
for iA = 1:size(aa,1)
if Sign == 1 && ~isempty(r{iA,1})% && a{iA,1}.IsPlural
continue; % Assume everything is not switchable
end
if Sign == -1 && isempty(r{iA,1})
continue;
end
for iS = 1:aa{iA,1}.nState
%%Membership function calculation probability
if iA == 7
p1 = 0;
end
p1 = Relia(AvgdP,aa{iA,1}.State{iS}.P,gStdP);%%haveproblem
p2 = Relia(AvgdQ,aa{iA,1}.Q,gStdP);
p3 = Relia(norm(AvgdIh),norm(aa{iA,1}.Ih),gStdI);
obj = + wP * p1 + wQ * p2 + wIh * p3;
if obj > selObj
selA = iA;
selS = iS;
selObj = obj;
selP = [p1 p2 p3];
end
end
end
if selObj < 0
disp('error happened');
%%Error logging
elseif selObj < 0.8 %%0.8 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%haveproblem20181006
fprintf('Result:Electrical identification complete, not electrombile.Info:%.3f£¬%.3f£¬%.3f',AvgdP,AvgdQ,AvgdIh);
elseif selObj <= 1
fprintf(' matching Type=%s, State=%s , its obj = %.3f , probability:%s \n',aa{selA,1}.Type , aa{selA,1}.State{selS,1}.StateName, selObj , selP);
sign = 1;
else
disp('error happened');
%%Error logging
end
if selA == 0
disp('Warning: There is something wrong !');
[PrePara,EventTkt] = ParaInit(PrePara,EventTkt);
continue;
end
if Sign == 1
SignStr = '[On]';
if selA > 0
r{selA,1} = 1;
end
else
SignStr = '[Off]';
if selA > 0
r{selA,1} = []; %%20181006
end
end
if selA > 0
fprintf(' - Event recogized: %s, %s, %s\n', SignStr, aa{selA,1}.Type, aa{selA,1}.State{selS}.StateName);
else
fprintf('No match data detected');
end
% wIh * norm(e{end,1}.dIHD-a{selA,1}.State{selS}.IHD));
% Update Status
% pause
figure(figAll);
LC_flag = 0;
if e{end,1}.dIr > 100
disp('Warning: Leakage Current is too large !');
subplot(4,3,figLeakCurrent); imshow(pic);drawnow;
LC_flag = 1;
elseif e{end,1}.dIr < -90
disp('Leakage Current is OK !');
subplot(4,3,figLeakCurrent); imshow(pic_gray);drawnow;
LC_flag = 0;
end
if strcmp(SignStr,'[On]')
switch aa{selA,1}.Type
case 'Microwave'
subplot(4,3,figWave) ;
case 'ECar'
subplot(4,3,figAuto) ;
case 'Refri'
subplot(4,3,figRefri) ;
case 'Cooker'
subplot(4,3,figCooker) ;
case 'Kettle'
subplot(4,3,figKettle) ;
case 'CFL'
subplot(4,3,figCFL) ;
case 'LED'
subplot(4,3,figLED) ;
case 'Hairdryer'
subplot(4,3,figHairdryer) ;
case 'Iron'
subplot(4,3,figIron) ;
case 'Heater'
subplot(4,3,figHeater) ;
case 'Bath'
subplot(4,3,figBath) ;
end
if ~LC_flag
imshow(pic_green);
else
imshow(pic);
end
elseif strcmp(SignStr , '[Off]')
switch aa{selA,1}.Type
case 'Wave'
subplot(4,3,figWave) ;
case 'Ecar'
subplot(4,3,figAuto) ;
case 'Refri'
subplot(4,3,figRefri) ;
case 'Cooker'
subplot(4,3,figCooker) ;
case 'Kettle'
subplot(4,3,figKettle) ;
case 'CFL'
subplot(4,3,figCFL) ;
case 'LED'
subplot(4,3,figLED) ;
case 'Hairdryer'
subplot(4,3,figHairdryer) ;
case 'Iron'
subplot(4,3,figIron) ;
case 'Heater'
subplot(4,3,figHeater) ;
case 'Bath'
subplot(4,3,figBath) ;
end
imshow(pic_gray);
end
drawnow;
end
end
end
end

Iniciar sesión para responder a esta pregunta.

Respuestas (1)

Mohammad Junayed
Mohammad Junayed el 25 de Oct. de 2018
here is the file

Categorías

Más información sobre Switches and Breakers en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by