Lstm 路面等级识别

路面等级分类

一般将路面等级分为A-H8个等级，车辆在不同路面等级下的路面响应不同，如车辆垂直加速度。

建立路面模型以及车辆1/4悬架模型,得到不同路面等级下的车辆簧下质量加速度输出。

数据集制作

随机截取某一路面等级下的加速度数据，记录该数据对应的路面等级,

clear
clc
for j=1:400
    if j<=100
    x=rand*20+8;
    Y_train(j,1)=1;
    elseif j>100&&j<=200
    x=rand*50+32;
    Y_train(j,1)=2;
    elseif j>200&&j<=300
    x=rand*300+128; 
    Y_train(j,1)=3;
    elseif j>300&&j<=400
    x=rand*1200+512;
    Y_train(j,1)=4;
    end
    y=fix(rand*20)+1;
    sqno_1=x*1e-6;                         %路面不平度系数
    u=2;                                   %车速
    noo=0.011;                             %空间截止频率
    f0=2*pi*noo*u;                         %下截止频率
    stepsize=1/105;
    time=1;
    simout=sim('xuanjiamox1_4_1',0:stepsize:time);  %运行simulink模型
    t=simout.t.data;
    gujizhi=simout.gujizhi_2.data;
    zt_dot2=simout.zt_dot2.data;
    zs_dot2=simout.zs_dot2.data;
    A1(j,1)=1;
    t=simout.t.data;
    gujizhi=simout.gujizhi_2.data;
    %% 利用burg算法求得路面位移功率谱
    [m,n]=size(gujizhi);
    p=fix(m/10);% p不存在，需要定阶，定阶准则即criterion,阶次一般不超过信号长度的1/3
    [psd, f, p] = myBurg(gujizhi,10, p);
    X_train{j,1}(1,:)=zt_dot2(y:y+50,1);
    X_train{j,1}(2,:)=zs_dot2(y:y+50,1); 
end
for j=1:400
    if j<=100
    x=rand*20+8;
    Y_test(j,1)=1;
    elseif j>100&&j<=200
    x=rand*50+32;
    Y_test(j,1)=2;
    elseif j>200&&j<=300
    x=rand*300+128; 
    Y_test(j,1)=3;
    elseif j>300&&j<=400
    x=rand*1200+512;
    Y_test(j,1)=4;
    end
    y=fix(rand*20)+1;
    sqno_1=x*1e-6;                         %路面不平度系数
    u=2;                                   %车速
    noo=0.011;                             %空间截止频率
    f0=2*pi*noo*u;                         %下截止频率
    stepsize=1/105;
    time=1;
    simout=sim('xuanjiamox1_4_1',0:stepsize:time);  %运行simulink模型
    t=simout.t.data;
    gujizhi=simout.gujizhi_2.data;
    zt_dot2=simout.zt_dot2.data;
    zs_dot2=simout.zs_dot2.data;
    A1(j,1)=1;
    t=simout.t.data;
    gujizhi=simout.gujizhi_2.data;
    %% 利用burg算法求得路面位移功率谱
    [m,n]=size(gujizhi);
    p=fix(m/10);% p不存在，需要定阶，定阶准则即criterion,阶次一般不超过信号长度的1/3
    [psd, f, p] = myBurg(gujizhi,10, p);
    X_test{j,1}(1,:)=zt_dot2(y:y+50,1);
    X_test{j,1}(2,:)=zs_dot2(y:y+50,1);
    %X_train{j,1}(3,:)=psd(y:y+50,1);    
end
save X_train.mat
save Y_train.mat
save X_test.mat
save Y_test.mat

X_train数据格式如下：n x 1 cell, 每个cell下包含2 x m 的数组，，2表示2种特征，分别为簧下加速度与簧上加速度,m表示数据长度，n表示样本数量

cell
2 x m double
2 x m double
2 x m double
...

Y_train数据格式为：n x 1

Lstm-路面等级识别实现

clear
clc
%% 数据加载
load X_train.mat
load Y_train.mat
load X_test.mat
load Y_test.mat
XTrain=X_train;
YTrain=categorical(Y_train);% double转换为categorical
%% LSTM设置
inputSize = 2;            %输入层数量
numHiddenUnits = 100;     %隐藏层数量
numClasses = 4;           %分类种类
layers = [ ...
    sequenceInputLayer(inputSize)
    bilstmLayer(numHiddenUnits,'OutputMode','last')
    fullyConnectedLayer(numClasses)
    softmaxLayer
    classificationLayer]
maxEpochs = 100;          %最大迭代次数
miniBatchSize = 100;
options = trainingOptions('adam', ...
    'ExecutionEnvironment','cpu', ...
    'GradientThreshold',4, ...
    'MaxEpochs',maxEpochs, ...
    'MiniBatchSize',miniBatchSize, ...
    'SequenceLength','longest', ...
    'Shuffle','never', ...
    'Verbose',0, ...
    'Plots','training-progress');
net = trainNetwork(XTrain,YTrain,layers,options);
%% 测试集测试
XTest=X_test;
YTest=categorical(Y_test);
YPred = classify(net,XTest, ...
    'MiniBatchSize',miniBatchSize, ...
    'SequenceLength','longest');
acc = sum(YPred == YTest)./numel(YTest)

训练结果：

测试集准确率：

序列至序列分类

随机生成不同等级路面数据

clear
clc
len1=30;
f = waitbar(0,'Please wait...');
for k=1:10
%% 训练集
    for j=1:len1
    nodes=[rand*200+300;rand*299+700];
    nodes=sortrows(nodes);
    for i=1:3
        num=rand*4 ;
        if num>0&&num<=1
           xx=rand*20+8;
           elseif num>1&&num<=2
           xx=rand*50+32;
           elseif num>2&&num<=3
           xx=rand*350+150;   
           elseif num>3&&num<=4
           xx=rand*1000+700;
           elseif num>4&&num<=5
           xx=rand*4500+3000; 
           elseif num>5&&num<=6
           xx=rand*20000+10000;
           elseif num>6&&num<=7
           xx=rand*100000+38000;
           elseif num>7&&num<=8
           xx=rand*300000+150000;
        end
        DJ(1,i)=xx*1e-6;
        dj(i,1)=fix(num);
    end
    DJ=DJ';
    u=2;                                   %车速
    noo=0.011;                             %空间截止频率
    f0=2*pi*noo*u;                         %下截止频率
    stepsize=1/10;
    time=1000;
    simout=sim('xuanjiamox1_4_1',0:stepsize:time);  %运行simulink模型
    Y_dj=simout.Y_dj.Data;
    zr=simout.zr.Data;
    zs_dot2=simout.zs_dot2.Data;
    zt_dot2=simout.zt_dot2.Data;
    zs=simout.zs.Data;
    zt=simout.zt.Data;
    zs_dot1=simout.zs_dot1.Data;
    zt_dot1=simout.zt_dot1.Data;
    X_train{k,1}{j,1}(1,:)=zr';
    X_train{k,1}{j,1}(2,:)=zs_dot2';
    X_train{k,1}{j,1}(3,:)=zt_dot2';
    X_train{k,1}{j,1}(4,:)=zs';
    X_train{k,1}{j,1}(5,:)=zt';
    X_train{k,1}{j,1}(6,:)=zs_dot1';
    X_train{k,1}{j,1}(7,:)=zt_dot1';
    Y_train{k,1}{j,1}=Y_dj';
    str=['训练集生成中...',num2str(100*j/len1),'%'];    
    waitbar(j/len1,f,str) ;
    end   
    %% 测试集
    len2=10;
    for j=1:len2
    nodes=[rand*200+300;rand*299+700];
    nodes=sortrows(nodes);
    for i=1:3
        num=rand*4 ;
        if num>0&&num<=1
           xx=rand*20+8;
           elseif num>1&&num<=2
           xx=rand*50+32;
           elseif num>2&&num<=3
           xx=rand*350+150;   
           elseif num>3&&num<=4
           xx=rand*1000+700;
           elseif num>4&&num<=5
           xx=rand*4500+3000; 
           elseif num>5&&num<=6
           xx=rand*20000+10000;
           elseif num>6&&num<=7
           xx=rand*100000+38000;
           elseif num>7&&num<=8
           xx=rand*300000+150000;
        end
        DJ(1,i)=xx*1e-6;
        dj(i,1)=fix(num);
    end
    DJ=DJ';
    u=2;                                   %车速
    noo=0.011;                             %空间截止频率
    f0=2*pi*noo*u;                         %下截止频率
    stepsize=1/10;
    time=1000;
    simout=sim('xuanjiamox1_4_1',0:stepsize:time);  %运行simulink模型
    Y_dj=simout.Y_dj.Data;
    zr=simout.zr.Data;
    zs_dot2=simout.zs_dot2.Data;
    zt_dot2=simout.zt_dot2.Data;
    zs=simout.zs.Data;
    zt=simout.zt.Data;
    zs_dot1=simout.zs_dot1.Data;
    zt_dot1=simout.zt_dot1.Data;
    X_test{k,1}{j,1}(1,:)=zr';
    X_test{k,1}{j,1}(2,:)=zs_dot2';
    X_test{k,1}{j,1}(3,:)=zt_dot2';
    X_test{k,1}{j,1}(4,:)=zs';
    X_test{k,1}{j,1}(5,:)=zt';
    X_test{k,1}{j,1}(6,:)=zs_dot1';
    X_test{k,1}{j,1}(7,:)=zt_dot1';
    Y_test{k,1}{j,1}=Y_dj';
    str=['测试集生成中...',num2str(100*j/len2),'%'];    
    waitbar(j/len2,f,str) ;
    end      
end
    save ('X_train.mat',"X_train")
    save ('Y_train.mat',"Y_train")
    save ('X_test.mat',"X_test")
    save ('Y_test.mat',"Y_test")

如图所示

LSTM对序列进行分类

clc
clear
load X_train
load Y_train
load X_test
load Y_test
for k=1:1%length(Y_train)
        len=length(Y_train{k});
        for i=1:len
        Y_train{k}{i} = categorical(Y_train{k}{i});
        end
        len=length(Y_test{k});
        for i=1:len
        Y_test{k}{i} = categorical(Y_test{k}{i});
        end
       %% 
        numFeatures = 7;
        numHiddenUnits = 50;
        numClasses = 4;
        miniBatchSize=128;
        layers = [ ...
            sequenceInputLayer(numFeatures)
            lstmLayer(numHiddenUnits,'OutputMode','sequence')
            fullyConnectedLayer(numClasses)
            softmaxLayer
            classificationLayer];
        options = trainingOptions('adam', ...
            'ExecutionEnvironment','gpu', ...
            'MaxEpochs',100, ...
            'MiniBatchSize',miniBatchSize,...
            'Shuffle','every-epoch',...
            'GradientThreshold',4, ...
            'InitialLearnRate',0.05,...
            'LearnRateSchedule','piecewise',...
            'LearnRateDropPeriod',30,...
            'Verbose',0, ...
            'Plots','training-progress');
        net = trainNetwork(X_train{k},Y_train{k},layers,options);
       %% test 
        for i=1:length(X_test)
        YPred{k}{1,i} = classify(net,X_test{k}{i});
        acc(1,i)= sum(YPred{k}{i} == Y_test{k}{i})./numel(Y_test{k}{i})
         figure
         confusionchart(Y_test{k}{i},YPred{k}{i})
        end
        Acc{k}=sum(acc)/length(X_test)
        %%
        for i=1:length(X_test{k})
        X=X_test{k}{i}(1,:);
        classes_test = categories(Y_test{k}{i});
        classes_pred = categories(YPred{k}{i});
        figure
            for j = 1:numel(classes_test)
                label_test = classes_test(j);
                label_pred = classes_pred(j);
                idx_test = find(Y_test{k}{i} == label_test);
                idx_pred = find(YPred{k}{i} == label_pred);
                hold on
                plot(idx_test,X(idx_test),".")
                plot(idx_pred,X(idx_pred),".")
            end
        end
end

结果如下：