%% ABOUT THIS FILE % ---------------------------------------------------------------------- % Type-Bp for KAGRA % Coded by A.Shoda on 2015/06/28 % ---------------------------------------------------------------------- %% PRELIMINARY clear all; % Clear workspace close all; % Close plot windows addpath('../../utility'); % Add path to utilities addpath('/usr/local/matlab/altmany-export_fig-76bd7fa'); g = 9.81; %% IMPORT SUSPENSION MODEL matfile='typeBp160620mdl'; load([matfile,'.mat']); %% TUNING DAMPER % This part compensates the failure in converting the structural damping to % viscous damping. % REDUCE DAMPING ON RIM sys1.a(34,34)=sys1.a(34,34)/30; % RIM sys1.a(40,40)=sys1.a(40,40)/30; % RRM sys1.a(46,46)=sys1.a(46,46)/30; % RTM % INCREASE DAMPING ON YF0 %sys1.a(15,15)=sys1.a(15,15)*3000; % YF1 %sys1.a(15,15)=sys1.a(15,15)*1000; % YF1 % INCREASE DAMPING ON LF1 %sys1.a(13,13)=sys1.a(13,13)*3000; % LF1 % INCREASE DAMPING ON TF1 %sys1.a(14,14)=sys1.a(14,14)*3000; % TF1 % INCREASE DAMPING ON YRM, YTM sys1.a(42,42)=sys1.a(42,42)*30; % YRM sys1.a(48,48)=sys1.a(48,48)*300; % YTM %% IMPORT SERVO FILTERS addpath('servofilter'); % Add path to servo typeBp160628_no_control; % NO CONTROL MODE rmpath ('servofilter'); % Remove path to servo %% IMPORT SIMULINK MODEL mdlfile='typeBpsimctrl_160620'; % typeB1 ver.150616 st =linmod(mdlfile); invl =strrep(st.InputName, [mdlfile,'/'],''); outvl =strrep(st.OutputName,[mdlfile,'/'],''); sysc0 =ss(st.a,st.b,st.c,st.d,'inputname',invl,'outputname',outvl); %% FREQUENCY freq1=logspace(-2,3,1001); freq2=logspace(-2,3,1001); freq=logspace(-2,2,1001); %% NOISE % NOISE data_OSEM = importdata('../../../../noise/OSEMnoiseworstproto_disp.dat'); % um/rtHz noise_OSEM = interp1(data_OSEM(:,1),data_OSEM(:,2)*1e-6,freq')'; %m/rtHz n_OSEM_LTM = noise_OSEM/sqrt(4); % [m/rtHz] n_OSEM_PTM = noise_OSEM/sqrt(2)/117e-3; % [rad/rtHz] n_OSEM_YTM = noise_OSEM/sqrt(2)/117e-3; % [rad/rtHz] n_OSEM_LIM = noise_OSEM/sqrt(1); % [m/rtHz] n_OSEM_TIM = noise_OSEM/sqrt(2); % [m/rtHz] n_OSEM_VIM = noise_OSEM/sqrt(3); % [m/rtHz] n_OSEM_RIM = noise_OSEM/sqrt(3/2)/96e-3;% [rad/rtHz] n_OSEM_PIM = noise_OSEM/sqrt(3/2)/96e-3;% [rad/rtHz] n_OSEM_YIM = noise_OSEM/sqrt(2)/87e-3; % [rad/rtHz] data_LVDT = importdata('../../../../noise/LVDTnoiseADC_disp.dat'); % um/rtHz noise_LVDT = interp1(data_LVDT(:,1),data_LVDT(:,2)*1e-6,freq')'; % m/rtHz %n_LVDT_LF0 = noise_LVDT/sqrt(3/2); % [m/rtHz] %n_LVDT_TF0 = noise_LVDT/sqrt(3/2); % [m/rtHz] %n_LVDT_YF0 = noise_LVDT/sqrt(3)/600e-3; % [rad/rtHz] %n_LVDT_VF0 = noise_LVDT*2; % [m/rtHz] n_LVDT_VF1 = noise_LVDT*2; % [m/rtHz] n_LVDT_VF2 = noise_LVDT*2; % [m/rtHz] n_LVDT_LF2 = noise_LVDT/sqrt(2); n_LVDT_TF2 = noise_LVDT/sqrt(2); n_LVDT_YF2 = noise_LVDT/sqrt(2)/300e-3; n_Oplev = zeros(size(freq))+1.e-7; %[rad/rtHz] (from Michimura's slide JGW-T1202403-v1) data_seism = importdata('../../../../noise/KamiokaSeismicHighNoise.dat'); % m/rtHz noise_seism = interp1(data_seism(:,1),data_seism(:,2),freq')';% m//rtHz %% NOISE MODEL % mypsdplotopt({noise_OSEM,noise_seism},freq,... % 'title','Noise Model',... % 'legend',{'OSEM','seismic'},... % 'color',{'r-','k-'},... % 'ylim',[1e-13,1e-5],'ylabel','Magnitude [m/rtHz] ') %% MODEL addpath('servofilter'); % Add path to servo typeBp160620_servo_damp; % Lock-acq MODE rmpath ('servofilter'); % Remove path to servo mdlfile='typeBpsimctrl_160620'; % typeBp ver.160620 st =linmod(mdlfile); invl =strrep(st.InputName, [mdlfile,'/'],''); outvl =strrep(st.OutputName,[mdlfile,'/'],''); sysc =ss(st.a,st.b,st.c,st.d,'inputname',invl,'outputname',outvl); %% CALC INITIAL [mag0Vg,~]=bodesus(sysc0,'accVGND','IFO_LTM',freq); [mag0S,~]=bodesus(sysc0,'accLGND','IFO_LTM',freq); mag0Vg=mag0Vg.*freq.*freq*pp*pp.*noise_seism; mag0S = mag0S.*freq.*freq*pp*pp.*noise_seism; mag0T=sumpsd({mag0Vg,mag0S}); % Mirror Displacement w/o Ctrl rms0T=makerms(freq,mag0T); % Mirror Displacement RMS w/o Ctrl rms0Tv=makerms(freq,mag0T.*freq*pp); % Mirror Velocity RMS w/o Ctrl %% CALC SEISM for LENGTH [magVg,~]=bodesus(sysc,'accVGND','IFO_LTM',freq); [magsS,~]=bodesus(sysc,'accLGND','IFO_LTM',freq); magVg=magVg.*freq.*freq*pp*pp.*noise_seism; magsS=magsS.*freq.*freq*pp*pp.*noise_seism; magSeisT=sumpsd({magsS,magVg}); % Seismic noise of mirror rmsSeisT=makerms(freq,magSeisT); % Seismic noise RMS of mirror rmsSeisTv=makerms(freq,magSeisT.*freq*pp); % Seismic noise velocity RMS of mirror %% CALC GAS for LENGTH %[magg1,~]=bodesus(sysc,'n_LVDT_VF1','IFO_LTM',freq); [magg2,~]=bodesus(sysc,'n_LVDT_VF2','IFO_LTM',freq); %magg1=magg1.*n_LVDT_VF1; magg2=magg2.*n_LVDT_VF2; magGasTot=sumpsd({magg2}); %GAS LVDT noise [magfL,~]=bodesus(sysc,'n_LVDT_LF2','IFO_LTM',freq); [magfT,~]=bodesus(sysc,'n_LVDT_TF2','IFO_LTM',freq); [magfY,~]=bodesus(sysc,'n_LVDT_YF2','IFO_LTM',freq); magfL=magfL.*n_LVDT_LF2; magfT=magfT.*n_LVDT_TF2; magfY=magfY.*n_LVDT_YF2; magLvdtTot=sumpsd({magfL,magfT,magfY}); %F2 LVDT noise %% CALC ABOUT OSEM for LENGTH [magoL,~]=bodesus(sysc,'n_OSEM_LIM','IFO_LTM',freq); [magoP,~]=bodesus(sysc,'n_OSEM_PIM','IFO_LTM',freq); [magoY,~]=bodesus(sysc,'n_OSEM_YIM','IFO_LTM',freq); [magoT,~]=bodesus(sysc,'n_OSEM_TIM','IFO_LTM',freq); [magoR,~]=bodesus(sysc,'n_OSEM_RIM','IFO_LTM',freq); [magoV,~]=bodesus(sysc,'n_OSEM_VIM','IFO_LTM',freq); magoL=magoL.*n_OSEM_LIM; magoT=magoT.*n_OSEM_TIM; magoV=magoV.*n_OSEM_VIM; magoR=magoR.*n_OSEM_RIM; magoP=magoP.*n_OSEM_PIM; magoY=magoY.*n_OSEM_YIM; [magtL,~]=bodesus(sysc,'n_OSEM_LTM','IFO_LTM',freq); [magtP,~]=bodesus(sysc,'n_OSEM_PTM','IFO_LTM',freq); [magtY,~]=bodesus(sysc,'n_OSEM_YTM','IFO_LTM',freq); magtL=magtL.*n_OSEM_LTM; magtP=magtP.*n_OSEM_PTM; magtY=magtY.*n_OSEM_YTM; magOsemTot=sumpsd({magoL,magoP,magoY,magoT,magoR,magoV,magtL,magtP,magtY}); %OSEM noise %% CALC ABOUT OPLEV for LENGTH [magolP,~]=bodesus(sysc,'n_OpLev_PTM','IFO_LTM',freq); [magolY,~]=bodesus(sysc,'n_OpLev_YTM','IFO_LTM',freq); magolP=magolP.*n_Oplev; %Assume Oplev noise level is same as OSEM magolY=magolY.*n_Oplev; magOlTot=sumpsd({magolP, magolY}); %% ALL Noise magnTot=sumpsd({magGasTot,magLvdtTot, magOsemTot, magOlTot}); %ALL Noise %magnTot=sumpsd({magGasTot, magOsemTot}); %ALL Noise rmsnTot=makerms(freq,magnTot); rmsnTotv=makerms(freq,magnTot.*freq*pp); %% ALL magAllT=sumpsd({magSeisT,magnTot}); % magaT=magoTot; rmsAllT=makerms(freq,magAllT); rmsAllTv=makerms(freq,magAllT.*freq*pp); %% Requirement reqPR=freq; for i=1:length(freq) if freq(i)<8 reqPR(i)=0; else reqPR(i)=5e-15*(freq(i)/10)^(-2)+5e-16; end end % dataBRSE=importdata('../../SuspensionRequirement/BRSE/DisplacementNoiseRequirementCopy5Hz.dat'); % dataBRSE=[dataBRSE(:,1),dataBRSE(:,2),dataBRSE(:,3),dataBRSE(:,4),dataBRSE(:,5),dataBRSE(:,6),... % dataBRSE(:,7),dataBRSE(:,8),dataBRSE(:,9)]; % reqPR_BRSE = interp1(dataBRSE(:,1),dataBRSE(:,8),freq')';% m//rtHz % % dataDRSE=importdata('../../SuspensionRequirement/DRSE/DisplacementNoiseRequirementCopy5Hz.dat'); % dataDRSE=[dataDRSE(:,1),dataDRSE(:,2),dataDRSE(:,3),dataDRSE(:,4),dataDRSE(:,5),dataDRSE(:,6),... % dataDRSE(:,7),dataDRSE(:,8),dataDRSE(:,9)]; % reqPR_DRSE = interp1(dataDRSE(:,1),dataDRSE(:,8),freq')';% m//rtHz %% PLOT % mypsdplotopt({magoTot,maggT,magsLG,magaT,mag0T,reqBS},freq,... % 'title','Control noise coupling during observation phase',... % 'legend',{'OSEM noise','GAS control noise','top control noise'... % 'With control total','No control',... % 'BS requirement'},... % 'color',{'m-','c-','g-','r-','b-','k--'},... % 'ylim',[1e-18,1e-5],'ylabel','Magnitude [m/rtHz] or [m/sec]') mypsdplotopt({magOsemTot,magGasTot,magLvdtTot,magOlTot,magAllT,mag0T,rmsAllT,rms0T,reqPR},freq,... 'title','Control noise coupling during damping phase',... 'legend',{'OSEM noise','GAS control noise','F2 LVDT control noise','Oplev control noise',... 'With control total','No control','With control total RMS','No control RMS','req. PR (BRSE)','req. PR (DRSE)'},... 'color',{'m-','g-','c-','y-','r-','b-','r--','b--','k-','k--'},... 'ylabel','Magnitude [m/rtHz]') % export_fig('figure/SensorNoiseCoupling/typeBp160620_Disp_damping.pdf') %% CALC INITIAL for Pitch [mag0pS,~]=bodesus(sysc0,'accLGND','PTM',freq); mag0pS=mag0pS.*freq.*freq*pp*pp.*noise_seism; % Pitch Angular Fluctuation w/o Ctrl rms0pT=makerms(freq,mag0pS); % Pitch RMS w/o Ctrl rms0pTv=makerms(freq,mag0pS.*freq*pp); %% CALC SEISM for Pitch [magpS,~]=bodesus(sysc,'accLGND','PTM',freq); magpS=magpS.*freq.*freq*pp*pp.*noise_seism; %% CALC LVDT for Pitch [magpfL,~]=bodesus(sysc,'n_LVDT_LF2','PTM',freq); [magpfT,~]=bodesus(sysc,'n_LVDT_TF2','PTM',freq); [magpfY,~]=bodesus(sysc,'n_LVDT_YF2','PTM',freq); magpfL=magpfL.*n_LVDT_LF2; magpfT=magpfT.*n_LVDT_TF2; magpfY=magpfY.*n_LVDT_YF2; magpLvdtTot=sumpsd({magpfL, magpfT, magpfY}); %% CALC OSEM for Pitch [magpoL,~]=bodesus(sysc,'n_OSEM_LIM','PTM',freq); [magpoP,~]=bodesus(sysc,'n_OSEM_PIM','PTM',freq); [magpoY,~]=bodesus(sysc,'n_OSEM_YIM','PTM',freq); [magpoT,~]=bodesus(sysc,'n_OSEM_TIM','PTM',freq); [magpoR,~]=bodesus(sysc,'n_OSEM_RIM','PTM',freq); [magpoV,~]=bodesus(sysc,'n_OSEM_VIM','PTM',freq); magpoL=magpoL.*n_OSEM_LIM; magpoP=magpoP.*n_OSEM_PIM; magpoY=magpoY.*n_OSEM_YIM; magpoT=magpoT.*n_OSEM_TIM; magpoR=magpoR.*n_OSEM_RIM; magpoV=magpoV.*n_OSEM_VIM; [magptL,~]=bodesus(sysc,'n_OSEM_LTM','PTM',freq); [magptP,~]=bodesus(sysc,'n_OSEM_PTM','PTM',freq); [magptY,~]=bodesus(sysc,'n_OSEM_YTM','PTM',freq); magptL=magptL.*n_OSEM_LTM; magptP=magptP.*n_OSEM_PTM; magptY=magptY.*n_OSEM_YTM; magpOsemTot=sumpsd({magpoL,magpoP,magpoY,magpoT,magpoR,magpoV,magptL,magptP,magptY}); %% CALC Oplev for Pitch [magpolP,~]=bodesus(sysc,'n_OpLev_PTM','PTM',freq); [magpolY,~]=bodesus(sysc,'n_OpLev_YTM','PTM',freq); magpolP=magpolP.*n_Oplev; %Assume that Oplev noise level is the same as OSEM magpolY=magpolY.*n_Oplev; magpOlTot=sumpsd({magpolP,magpolY}); %% ALL magpAllT=sumpsd({magpS,magpLvdtTot,magpOsemTot,magpOlTot}); rmspAllT =makerms(freq, magpAllT); rmspAllTv=makerms(freq, magpAllT.*freq*pp); %% PLOT PITCH mypsdplotopt({mag0pS,magpLvdtTot,magpOsemTot,magpOlTot,magpAllT,rms0pT,rmspAllT},freq,... 'title','Angular Fluctuation of the Mirror',... 'legend',{'No control','F2 LVDT noise',... 'OSEM noise','Oplev noise','total with control','no control RMS',... 'with control RMS'},... 'color',{'b-','k-','g-','y-','r-','c--','m--'},... 'ylim',[1e-12,1e-4],'ylabel','Magnitude [rad/rtHz] or [rad]') export_fig('figure/SensorNoiseCoupling/typeBp160620_Pitch_damping.pdf') %% PLOT Velocity magAllTv = magAllT.*freq*pp; mag0Tv = mag0T.*freq*pp; %rms0Tv = makerms(freq, mag0Tv); mypsdplotopt({mag0Tv,magOsemTot.*freq*pp,magGasTot.*freq*pp,magLvdtTot.*freq*pp,magOlTot.*freq*pp,magAllTv,rms0Tv,rmsAllTv},freq,... 'title','Velocity of the Mirror',... 'legend',{'No control','OSEM noise', 'GAS noise', 'F2 LVDT noise','Oplev noise',... 'With control total',... 'No control RMS','With control RMS'},... 'color',{'b-','k-','k--','g-','y-','r-','c--','m--'},... 'ylim',[1e-12,1e-4],'ylabel','Magnitude [m/s/rtHz] or [m/s]') % export_fig('figure/SensorNoiseCoupling/typeBp160620_Vel_damping.pdf') %% % [mag0pS,~]=bodesus(sysc0,'accLGND','OpLev_PTM',freq); % mag0pS=mag0pS.*freq.*freq*pp*pp.*noise_seism; % rms0pT=makerms(freq,mag0pS); % mypsdplotopt({mag0pS,rms0pT},freq,... % 'title','Noise coupling from OSEMs',... % 'legend',{... % 'Suspended mirror w/o Ctrl','Suspended mirror RMS'},... % 'color',{'r-','m--'},... % 'ylim',[1e-12,1e-4],'ylabel','Angular Fluctuation [rad/rtHz]') %% Force on F2 [magLf2,~]=bodesus(sysc,'accLGND','act_LF2_mon',freq); magLf2 = magLf2.*freq.*freq.*pp.*pp.*noise_seism; rmsLf2 = makerms(freq, magLf2); mypsdplotopt({magLf2, rmsLf2},freq, 'title', 'Feedback Force on BF',... 'legend',{'Feedback force', 'RMS'}, 'color',{'b','c'},... 'ylim', [1e-10,1e-2],'ylabel','Force [N/rtHz] or [N]'); export_fig('figure/SensorNoiseCoupling/typeBp160620_ForceOnF2_damping.pdf'); %%Plot 1/e damping time with damping controls ON amd OFF getqt_damp(sysc0,sysc) % export_fig('figure/SensorNoiseCoupling/typeBp160620_Qfactor_damping.pdf') %saving point on 2016.6.21