function [fm,mg,f,OLG]=CompareMeasGain(filename);
%Load measurement and compare with model

S1darm
S1darm_08
[A,B,C,D] = linmod('S1darm_08');
darmsys = ss(A,B,C,D);

%filename = '020824_0818';

eval(['load ' filename 'coh.txt']);
eval(['load ' filename 'tf.txt']);

eval(['Xcoh=X' filename 'coh;']);
eval(['Xtf=X' filename 'tf;']);

fm = Xcoh(:,1);
coh= Xcoh(:,2);
mg = (Xtf(:,2)+i*Xtf(:,3)); % meas. (ASQ/EXC)/(DARM_ERR/EXC)
mg = -asq2darm*mg;           %measured open loop gain
cg=1./(1+mg);               %measured closed loop gain

fg = (Xtf(:,4)+i*Xtf(:,5)); % meas. (DARM_CTRL/EXC)/(DARM_ERR/EXC)

gc=find(coh>.5); %good coherence points
f1=fm(min(gc));
f2=fm(max(gc));

f=logspace(log10(f1),log10(f2),1e4);
F=mybodesys(darmsys(2,2),f);
CLG=mybodesys(darmsys(1,1),f);
OLG=-F./CLG;
ugf=f(min(find(abs(OLG)<1)));  
%modelled open loop gain 
OLGm=interp1(f,OLG,fm);     
%filter gain from DARM_ERR to DARM_CTRL
FG=mybodesys(DARMDF*GW_K,fm);

figure(1)
clf

subplot(221)
semilogx(f,20*log10(abs(OLG)),fm(gc),20*log10(abs(mg(gc))),'r*', ...
	 'Linewidth',2)

i1=min(find(f>=f1));
i2=min(find(f>=min(f2,3e3)));


axis([f1 min(f2,3e3) -40 20])
ylabel('mag (dB)')
text(100,25,['Open Loop Gain ' filename ', ugf=' num2str(ugf,4) ...
       ' Hz, al\_gain= ' num2str(al_gain,3)])
legend('modelled','measured',3)
grid on

subplot(223)
semilogx(f,angle(OLG)*180/pi,fm(gc),angle(mg(gc))*180/pi,'r*','Linewidth',2)
axis([f1 min(f2,3e3) -180 180])
set(gca,'YTick',[-180:90:180])
ylabel('phase (deg)')
xlabel('Hertz')
grid on

subplot(222)
semilogx(fm(gc),abs(OLGm(gc)./mg(gc)),'.-')
axis([f1 min(f2,3e3) .5 1.5])
set(gca,'YTick',[.6:.2:1.4])
ylabel('fit/meas (abs)')
grid on


subplot(224)
semilogx(fm(gc),(180/pi)*angle(OLGm(gc)./mg(gc)),'.-')
axis([f1 min(f2,3e3) -10 10])
ylabel('fit/meas (phase)')
grid on
shg


figure(2)
clf


subplot(211)
semilogx(fm,20*log10(abs(FG)),fm,20*log10(abs(fg)),'r*',...
	 'Linewidth',2)
ax=axis;
axis([f1 min(f2,3e3) ax(3:4)])
ylabel('mag (dB)')
title(['DARM\_CTRL/DARM\_ERR ' filename])
legend('modelled','measured',4)
grid on

subplot(212)
semilogx(fm,angle(FG)*180/pi,fm,angle(fg)*180/pi,'r*','Linewidth',2)
axis([f1 min(f2,3e3) -180 180])
set(gca,'YTick',[-180:90:180])
ylabel('phase (deg)')
xlabel('Hertz')
grid on
shg

figure(3)
clf



subplot(211)
semilogx(f,20*log10(abs(CLG)),fm(gc),20*log10(abs(cg(gc))),'r*', ...
	 'Linewidth',2)
grid on

i1=min(find(f>=f1));
i2=min(find(f>=min(f2,3e3)));


axis([f1 min(f2,3e3) -20 10])
ylabel('mag (dB)')
title(['Closed Loop Gain ' filename])
legend(['modelled al\_gain= ' num2str(al_gain,3)],'measured',4)


subplot(212)
semilogx(f,angle(CLG)*180/pi,fm(gc),angle(cg(gc))*180/pi,'r*','Linewidth',2)
axis([f1 min(f2,3e3) -180 180])
set(gca,'YTick',[-180:90:180])
ylabel('phase (deg)')
xlabel('Hertz')
grid on
shg


keyboard

%Rana's entry, elog LLO 8/9/02
%Results (in nm/count @ DC):
%
%ITMX = 2.78
%ITMY = 2.88
%ETMX = 2.53
%ETMY = 2.67
%
%CARM = DARM = 1.72


