Here's .mat file illustrating above question.

Thank you,
Ted

• https://files.engineering.com/getfile.aspx?folder=6f250554-8b7d-48fa-8723-7

Here's live script file

Thank you,
Ted

• https://files.engineering.com/getfile.aspx?folder=0c1e6bde-f717-44ca-ae20-3

Quote:

decibels to magnitude

A decibel is a ratio; you need to find the reference level to convert to a magnitude

TTFN (ta ta for now)
I can do absolutely anything. I'm an expert! https://www.youtube.com/watch?v=BKorP55Aqvg
FAQ731-376: Eng-Tips.com Forum Policies forum1529: Translation Assistance for Engineers Entire Forum list http://www.eng-tips.com/forumlist.cfm

Thank you, IRstuff.
Would you mind if I ask you to provide a little more detail regarding finding reference level?

Thank you,
Ted

That's a question for the shaker lab.

TTFN (ta ta for now)
I can do absolutely anything. I'm an expert! https://www.youtube.com/watch?v=BKorP55Aqvg
FAQ731-376: Eng-Tips.com Forum Policies forum1529: Translation Assistance for Engineers Entire Forum list http://www.eng-tips.com/forumlist.cfm

Hmm...
Let me try again. I think that my poor English didn't allow me to express my question properly.
As far as I understand (please correct me if I'm wrong). Db is a ratio of magnitude of two signals expressed as 20*Log(A2/A1). Is there a way that I can get the results of tfestimate function to be expressed on linear scale?

Thank you,
Ted

Sure... Just undo the logarithm

TTFN (ta ta for now)
I can do absolutely anything. I'm an expert! https://www.youtube.com/watch?v=BKorP55Aqvg
FAQ731-376: Eng-Tips.com Forum Policies forum1529: Translation Assistance for Engineers Entire Forum list http://www.eng-tips.com/forumlist.cfm

Thank you, IRstuff,
Problem is that this is my 7th day with MatLab and I know very little at this point. It seems to me that conversion to Db is sort of imbedded in tfestimate MatLab function and I can't figure out how to undo the logarithm in the result plot of this function. When I write the function it may look like this: tfestimate(y,x,[],[],[],Fs) in brackets I can input window size, overlap, etc. But I don't see Log10 in input arguments options and so I can't choose not to enter it... I suspect that I'm missing something very simple but I can't find it in help.
May I ask you to suggest how exactly I can write it in MatLab to get rid of Logarithm?

Thank you,
Ted

The dB reference value is 1, so if you were testing a spring mass system with a k of 100 it would start at 20 dB at 0 Hz.

It is a good idea to do the fft of each channel separately and then divide one by the other to check exactly what is going on scaling wise. Alternatively generate two known signals (say A and 10*A) and see what tfestimate gives.

Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376: Eng-Tips.com Forum Policies http://eng-tips.com/market.cfm?

Having said that read the examples page for tfestimate, it does not automatically plot dB. https://www.mathworks.com/help/signal/ref/tfestima...

Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376: Eng-Tips.com Forum Policies http://eng-tips.com/market.cfm?

Quote:

The dB reference value is 1

It's "1" of something, but the OP should already know that, since the transfer function required dividing the response by the input, so all the OP needs to answer is WHAT is the input?

TTFN (ta ta for now)
I can do absolutely anything. I'm an expert! https://www.youtube.com/watch?v=BKorP55Aqvg
FAQ731-376: Eng-Tips.com Forum Policies forum1529: Translation Assistance for Engineers Entire Forum list http://www.eng-tips.com/forumlist.cfm

Thank you Greg,
There's definately a scaling mess which I have to fix first.
Regarding your 2nd reply - do you mean plotting abs of tfestimate? I didn't find any other way. May be it's my English playing tricks with me but I couldn't find any other way of plotting other than Db on y axis.

Thank you,
Ted

Yes, the default is to display dB, but if you take the output of tfestimate you can display that vector however you like, it's just a vector of complex numbers. ie

tfout=tfestimate(usual stuff in here);
plot(real(tfout),imag(tfout))

will give you a nyquist plot

Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376: Eng-Tips.com Forum Policies http://eng-tips.com/market.cfm?

Thank you Greg,
That's what I'm doing now and with some fiddling with signal(exported from Adams) scaling things are starting to make sense.

Thank you,
Ted

In fact I was doing it differently:
tfout=tfestimate(usual stiff in here)
tfMag=abs(tfout)
tfPhase=angle(tfout)*57.3 to approximate phase in deg.
Then plot the above as two separate plots.
After rescaling Adams signals this gives me reasonable results with heave.
Doing it the way you suggested returns a very strange plot (sorry).
I will see how things work with other signals (pitch and cpl in particular) and report.

Thank you,
Ted

BTW: how is your porpoising problem study going ?

So: here's what I do for computing lateral handling analysis:
I have data channels from road & sim data with the usual names. They are segments of length power of 2. Beginning and ending with a zero input for a second or 2.
nfft = pow2(round(log2(length(STEER))));
sps = 100; % sample rate
nsegs = size(STEER,2);
for nseg=1 %:nsegs
[sstxy,f] = tfestimate(STEER(:,nseg), LATACC(:,nseg),rectwin(nfft),nfft/2,nfft,sps); % ayg by steer transfer function
[rolltxy,f] = tfestimate(LATACC(:,nseg),ROLL (:,nseg),rectwin(nfft),nfft/2,nfft,sps); % roll by ay transfer function
[rollms2txy,f] = tfestimate(LATACC(:,nseg)*9.806,ROLL (:,nseg),rectwin(nfft),nfft/2,nfft,sps); % roll by ay transfer function
[betatxy,f] = tfestimate(LATACC(:,nseg),SIDSLP(:,nseg),rectwin(nfft),nfft/2,nfft,sps); % beta by ay transfer function
[betasteertxy,f]= tfestimate(STEER (:,nseg),SIDSLP(:,nseg),rectwin(nfft),nfft/2,nfft,sps); % beta by steer transfer function
[yawvtxy,f] = tfestimate(STEER(:,nseg), YAWVEL(:,nseg),rectwin(nfft),nfft/2,nfft,sps); % yawv by steer transfer function
[yawatxy,f] = tfestimate(STEER(:,nseg), YAWACC(:,nseg),rectwin(nfft),nfft/2,nfft,sps); % yaw accel by steer transfer function
[rollyawvtxy,f] = tfestimate(YAWVEL(:,nseg),ROLL (:,nseg),rectwin(nfft),nfft/2,nfft,sps); % roll by yawvel transfer function
[rollswatxy,f] = tfestimate(STEER(:,nseg) ,ROLL (:,nseg),rectwin(nfft),nfft/2,nfft,sps); % roll by steer transfer function

inx=find(f > 3.0 | f < .02);
f(inx)=[];
sstxy(inx)=[];
rolltxy(inx)=[];
rollms2txy(inx)=[];
betatxy(inx)=[];
betasteertxy(inx)=[];
yawvtxy(inx)=[];
yawatxy(inx)=[];
rollyawvtxy(inx)=[];
rollswatxy(inx)=[];

delta_f=f(2)-f(1);

handles.ss_gain = 100.*abs(sstxy); % convert to Bode form for gain
handles.ss_phase = angle(sstxy).*180/pi;
handles.yawv_gain = 100.*abs(yawvtxy);
handles.yawv_phase = angle(yawvtxy).*180/pi;
handles.yawa_gain = 100.*abs(yawatxy);
handles.yawa_phase = angle(yawatxy).*180/pi;
handles.roll_gain = abs(rolltxy);
handles.rollms2_gain = abs(rollms2txy);
handles.roll_phase = angle(rolltxy).*180/pi;
handles.rollms2_phase = angle(rollms2txy).*180/pi;
handles.beta_gain = abs(betatxy);
handles.beta_phase = angle(betatxy).*180/pi;
handles.betaswa_gain = 100*abs(betasteertxy);
handles.betaswa_phase = angle(betasteertxy).*180/pi;
handles.rollswa_gain = 100*abs(rollswatxy);
handles.rollswa_phase = angle(rollswatxy).*180/pi;
handles.rollyawv_gain = abs(rollyawvtxy);
handles.rollyawv_phase = angle(rollyawvtxy).*180/pi;
handles.f = f;
end

Thanks for great post Cibachrome,
I’m traveling right now but will reply as soon as I’m back home.

Thank you,
Ted

%Try this for the heave response with your data: Seems to work pretty well
load Question
[heavetxy,freq] = tfestimate(act_dz,heave,rectwin(nfft),nfft/2 ,nfft,Fs); % heave by act_dz transfer function
inx=find(freq > 20 | freq < .01);
freq(inx)=[];
heavetxy(inx)=[];
heave_gain = abs(heavetxy);
heave_phase = angle(heavetxy).*180/pi;

figure('name','T_M_S Porpoising Study','numbertitle','off','menubar','none')
subplot(2,2,1)
plot(freq,heave_gain,'o')
ylabel('Heave Gain (units ?)')
grid on
title('This is all you need to do Ted:')
subplot(2,2,3)
plot(freq,heave_phase,'o')
grid on
xlabel('Frequency (Hz.)')
ylabel('Phase (deg)')
sidetext('BillCobb zzvyb6@yahoo.com')

% Now go get the transfer function:
options = optimset('MaxFunEvals',1000000,'TolFun',.0000001,'Display','on');
FM0=[ .1 1 0.15 0.23 ]
for n=1: 5
[FM,error] = lsqcurvefit('splane2_fit',FM0 ,freq ,heave_gain,[],[],options);
heavesys =tf([FM(1:2 )],[FM(3:4) 1 ])
FM0=FM+ .001*rand(1);
end

N_1=num2str(FM(1));
N_0=num2str(FM(2));
D_2=num2str(FM(3));
D_1=num2str(FM(4));
D_0= '1';

f=.0 :.02:20;
[mag,phase]=bode(heavesys,f);
subplot(2,2,1)
hold on
plot(f,squeeze(mag),'b-')
mytexstr = ['$\frac{H(s)}{Z(s)}$ = $\frac{' N_1 's+' N_0 '}{' D_2 's^2+' D_1 's+' D_0 '}$']
text('string',mytexstr,'interpreter','latex','fontsize',15,'units','norm','pos',[.25 .75],'Fontweight','Bold');

[wn,z]=damp(heavesys);
subplot(2,2,3)
hold on
plot(f,squeeze(phase),'b-')
mytexstr = [ '$\omega_n = ' num2str(wn(1),3) ' Hz $']
text('string',mytexstr,'interpreter','latex','fontsize',12,'units','norm','pos',[.4 .7],'Fontweight','Bold');
mytexstr = [ '$\zeta = ' num2str(z(1),3) '$']
text('string',mytexstr,'interpreter','latex','fontsize',12,'units','norm','pos',[.4 .5],'Fontweight','Bold');

subplot(2,2,2)
bode(heavesys)
ylabel('Heave Response')
grid on

subplot(2,2,4)
step(heavesys,5)
ylabel('Heave Response (units ??)')
grid on

You'll need this: Note: I'm only using magnitude for fit criteria. phase floats
function [mag] = splane2_fit(FM,f)
sys = tf([FM(1:2 )],[FM(3:4) 1]);
mag = squeeze(bode(sys,f));

• https://files.engineering.com/getfile.aspx?folder=20599d90-0060-47cf-bc2e-9

Anybody considered that the porpoising is a ride problem, not necessarily an aero problem ? Sure, aero can attenuate it, but the root cause seems to be something else.

BTW: Your 1000 Hz. sample rate is way overkill. What would help your activity is a longer data segment for a smaller delta_F frequency increment, multiple segments for ensemble averaging, adding a starting & ending signal of a second or 2, and a power of 2 number of samples. Some nonlinear effects smear the FFT data from the Question.mat file, but hopefully, you get the idea...

This is pitch by act_dz.

• https://files.engineering.com/getfile.aspx?folder=55c9e819-1bba-4dca-a4b4-8

In Newey's book he mentions it as being a problem with an early Leyton House car. Basically as the front tire flexed the front of the floor got too close the ground , blocking the flow, or stalling the airflow, so it lost downforce, the front came up, and then the whole cycle repeated itself. Hard to fix with a passive suspension because the main compliance is undamped.

Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376: Eng-Tips.com Forum Policies http://eng-tips.com/market.cfm?

OK, I'm back home and able to reply.
Many thanks for your help Cibachrome!
I already found answers to my initial questions - it was a matter of getting the units/scaling straight:

For a frequency response I was using (Live script):

Heave_ Chassis Acceleration Response Magnitude Ratio

[h_az_tf,n]=tfestimate(act_az,heave_az,winsize,overlap,nfft,Fs);
[h_az_tf2,n]=tfestimate(act_az2,heave_az2,winsize,overlap,nfft,Fs);
plot(n,abs(h_az_tf),'r')
xlim(l)
hold on
plot(n,abs(h_az_tf2),'b')
hold off
title("Heave Chassis Acceleration FRF Ratio")
xlabel('Frequency (Hz)')
ylabel('Magnitude')
legend('d1',"d2")

Heave_ Chassis Acceleration Response Phase

plot(n,angle(h_az_tf)*57.3,'r')
xlim(l)
ylim([-200 200])
hold on
plot(n,angle(h_az_tf2)*57.3,'b')
hold off
title("Heave Chassis Displacement FRF PHASE")
xlabel('Frequency (Hz)')
ylabel('Phase Angle')
legend('d1',"d2")

d1_HeaveAcc_FRF_Max=max(abs(h_az_tf))
d2_HeaveAcc_FRF_Max=max(abs(h_az_tf2))
d1_HeaveAcc_FRF_RMS=rms(abs(h_az_tf))
d2_HeaveAcc_FRF_RMS=rms(abs(h_az_tf2))


This was to compare two different damper settings.
I fugured out that I need to increse the simulation time to get better resolution at each frequency step - this seems to be inline with your suggestion from your latest post.

Now let me study your code a little. This is a bit above my current MatLab knowledge level (remember my MatLab experience is less than a month at this point).

Thank you,
Ted

PS: Attached is Acar simulation data files

• https://files.engineering.com/getfile.aspx?folder=8fd4e693-acf2-48b7-b089-e

And here's live script files.

• https://files.engineering.com/getfile.aspx?folder=2a66db9e-54b7-4e22-9d45-9

Here's processing live script. There are also similar live scripts for Pitch, CPLV and tire deflection.

Thank you,
Ted

• https://files.engineering.com/getfile.aspx?folder=89060b77-7edd-4c02-a54d-d

There is a Matlab Controls setup panel that is supposed to allow you to choose units and display schemes.

ctrlpref

But it doesn't seem to be working for me. Set the preferences and reload Matlab.

But, the danger is ruining all your previous codes ! (as it did mine !)

Thanks Cibachrome,
Thing is that requests in A/Car are predefined and it is a bit too complicated to redefine them atm. So I set the units/scaling manually as can be seen in DataImport.mlx that I uploaded.


Thank you,
Ted

So far I'm having problem with adopting your transfer function code to my current data (see above uploaded .mat file).
This portion:

% Now go get the transfer function:
options = optimset('MaxFunEvals',1000000,'TolFun',.0000001,'Display','on');
FM0=[ .1 1 0.15 0.23 ]
for n=1: 5
[FM,error] = lsqcurvefit('splane2_fit',FM0 ,freq ,heave_gain,[],[],options);
heavesys =tf([FM(1:2 )],[FM(3:4) 1 ])
FM0=FM+ .001*rand(1);
end

I added above code to this live script:

[h_az_tf,n]=tfestimate(act_az,heave_az,winsize,overlap,nfft,Fs);
[h_az_tf2,n]=tfestimate(act_az2,heave_az2,winsize,overlap,nfft,Fs);
plot(n,abs(h_az_tf),'r')
xlim(l)
hold on
plot(n,abs(h_az_tf2),'b')
hold off
title("Heave Chassis Acceleration FRF Ratio")
xlabel('Frequency (Hz)')
ylabel('Magnitude')
legend('d1',"d2")

plot(n,angle(h_az_tf)*57.3,'r')
xlim(l)
ylim([-200 200])
hold on
plot(n,angle(h_az_tf2)*57.3,'b')
hold off
title("Heave Chassis Displacement FRF PHASE")
xlabel('Frequency (Hz)')
ylabel('Phase Angle')
legend('d1',"d2")
heave_gain=abs(h_az_tf)

Here I add your code with changes apropriate for my code above

options = optimset('MaxFunEvals',1000000,'TolFun',.0000001,'Display','on');
FM0=[ .1 1 0.15 0.23 ]
for n=1: 5
[FM,error] = lsqcurvefit('splane2_fit',FM0 ,n ,abs(h_az_tf),[],[],options);
heavesys =tf([FM(1:2 )],[FM(3:4) 1 ])
FM0=FM+ .001*rand(1);
end

N_1=num2str(FM(1));
N_0=num2str(FM(2));
D_2=num2str(FM(3));
D_1=num2str(FM(4));
D_0= '1';

Problem that it gives me an error:
Error using lsqcurvefit (line 286)
Function value and YDATA sizes are not equal.

How do I get your function splane2_fit to be the same size as my YDATA (abs(h_az_tf)?

Thank you,
Ted

"Thing is that requests in A/Car are predefined and it is a bit too complicated to redefine them atm"

I nearly fell of my chair at that one. In order to make our custom scripts work with /Car, MSC (or some lucky contractor) have written a black box that takes Car req files and re-mangles them back into the form we used to have in Chassis, where REQ numbers were fixed, for example REQ/1080 was always toe caster camber. I am not aware of the gory details.

Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376: Eng-Tips.com Forum Policies http://eng-tips.com/market.cfm?