200 lines
		
	
	
		
			5.2 KiB
		
	
	
	
		
			Matlab
		
	
	
	
	
	
			
		
		
	
	
			200 lines
		
	
	
		
			5.2 KiB
		
	
	
	
		
			Matlab
		
	
	
	
	
	
| %% Title: id31 microstation in EXP hutch for RAW data
 | |
| % Date: 12 october 2018
 | |
| 
 | |
| %% NOTE
 | |
| % With this file you can load the raw mat files
 | |
|  
 | |
| %% Description: measure on id31 microstation in exp hutch
 | |
|  
 | |
| % FS: =256Hz
 | |
| 
 | |
| 
 | |
| 
 | |
| % all L28 (31V/m/s) except CH1 L4-C (276V/m/s)
 | |
| % ch1: marble Z
 | |
| % ch2: outer frame Ty Z
 | |
| % ch3: inner frame Tilt Z
 | |
| % ch4: hexa Z
 | |
| % ch5: marble H
 | |
| % ch6: outer frame TY H
 | |
| % ch7: inner frame Tilt H
 | |
| % ch8: hexa H
 | |
| % ch9: hammer
 | |
| 
 | |
| %% measurements 12 october 2018
 | |
| 
 | |
| % excitation Y marble corner
 | |
| % Measurement1
 | |
| 
 | |
| % excitation Y outer frame corner
 | |
| % Measurement2
 | |
| 
 | |
| % excitation Y hexa
 | |
| % Measurement3
 | |
| 
 | |
| % ------------------------------------------------
 | |
| 
 | |
| % excitation Z marble corner
 | |
| % Measurement4
 | |
| 
 | |
| % excitation Z outer frame corner
 | |
| % Measurement5
 | |
| 
 | |
| % excitation Z hexa
 | |
| % Measurement6
 | |
| 
 | |
| % ------------------------------------------------
 | |
| 
 | |
| % excitation X marble corner
 | |
| % Measurement7
 | |
| 
 | |
| % excitation X outer frame corner
 | |
| % Measurement8
 | |
| 
 | |
| % excitation X hexa
 | |
| % Measurement9
 | |
| 
 | |
| %%
 | |
| microstation=['Marble     '; 'TY         ';'Tilt       '; 'Hexapod    '];
 | |
| 
 | |
| % 
 | |
| %% PARAMETERS 
 | |
| 
 | |
| beamline='ID31 Nanostation - Hammer testing';
 | |
| % --------------------------------
 | |
| %%----------OROS -----------------
 | |
| ch_max=16;
 | |
| % --------------------------------
 | |
| 
 | |
| %mult=1e6/276*173;  % --> m/s to micron/s and sensitivity correction for L4-C
 | |
| 
 | |
| nyqhp=2.56; % nyquist 
 | |
| f_cut=0.5; % cut frequency for high pass filter
 | |
| t_win=4; % window length in sec
 | |
| t_ovlp=0; % overlap window in sec
 | |
| 
 | |
| %d=1.3; % distance between vertical sensors
 | |
| 
 | |
| 
 | |
| warning off MATLAB:divideByZero
 | |
| 
 | |
| % specify capt # for which to run this
 | |
| capt=[1:9];
 | |
| 
 | |
| % specify channels for which shut correction must be applied
 | |
| shunt_ch=[1];
 | |
| 
 | |
| % in case of hammer inpacts specify capt # where it doesnt occur
 | |
| no_hammer=[];
 | |
| %no_hammer=[0];
 | |
| % specify hammer channel (or ch to find peak due to impacts)
 | |
| shock_ch=9;
 | |
| 
 | |
| %% main loop --------
 | |
| % ------------------
 | |
| for i=capt
 | |
|     
 | |
|         
 | |
|     eval(['load Measurement_raw',num2str(i)])
 | |
|     freq_max=Track1_TrueBandWidth;
 | |
|     dts=1/(freq_max*nyqhp); 
 | |
|     
 | |
|     freq=linspace(0,freq_max,t_win*freq_max);
 | |
|     wo=2*pi*freq;
 | |
|     
 | |
|     for k=1:ch_max
 | |
|     vname=['Track',num2str(k)];
 | |
|     array_exist(k)=ismember(vname,who);
 | |
|     end
 | |
|     non_zero=find(array_exist);
 | |
|     for z=non_zero(1):length(non_zero)
 | |
|     track_nb=['Track',num2str(z)]'; 
 | |
|     eval(['data(:,z)=Track',num2str(z),';']); 
 | |
|     end
 | |
|     c=data*mult;
 | |
|     
 | |
| 
 | |
|         
 | |
|     %-------------
 | |
|     nbch=size(c,2);
 | |
|     %-------------
 | |
|     r=length(c);
 | |
|     if r/2~=fix(r/2)    % loop to test for odd or even nb of samples
 | |
|        c=c(1:r-1,:);   % take only even
 | |
|     else 
 | |
|     end
 | |
|     %------------------------------
 | |
|     time=linspace(0,length(c)*dts,length(c)); 
 | |
|     
 | |
|     for j=nbch %shunt_ch
 | |
|         [c(:,j),c_shut]=shut_c(c(:,j),1/dts);   % correct for shunt 
 | |
|     end
 | |
|     %c(:,8)=(c(:,5)-c(:,4))/d; % differential Theta Y angle
 | |
|     
 | |
|      b=find(no_hammer==i); %      if i==1 | i==2 | i==6 
 | |
|      if b~=0
 | |
|          [psd_v,integ_v,psd_d,integ_d]=integrated_psd(c,t_win,t_ovlp,nyqhp,dts);
 | |
|          [frz_cut,crsp,pwsp,coherz,nsp]=fqresp(c,shock_ch,t_win,t_ovlp,nyqhp,dts);
 | |
|          
 | |
|      else 
 | |
|          
 | |
|         thresh=0.5;  % threshold of max value
 | |
|         sep=2.5;  % separation minimum of peaks in sec
 | |
|         pre_ev=2; % pre event delay in sec
 | |
|         pos_ev=2; % post event delay in sec
 | |
|                 
 | |
|         [ti,t_impact]=findpeaks(c(:,shock_ch),'minpeakheight',max(c(:,shock_ch))*thresh,'minpeakdistance',ceil(sep/dts));
 | |
|         % find times at which there are impacts (threshold of max and separated by sep sec)
 | |
|         
 | |
|         psd_v=zeros((pre_ev+pos_ev)/dts/nyqhp,nbch);
 | |
|         psd_d=zeros((pre_ev+pos_ev)/dts/nyqhp,nbch);
 | |
|         frz_cut=zeros((pre_ev+pos_ev)/dts/nyqhp,nbch);
 | |
|        
 | |
|         for k=1:length(t_impact)
 | |
|             ibeg=fix(t_impact(k)-(pre_ev/dts));
 | |
|             iend=fix(t_impact(k)+(pos_ev/dts));
 | |
|             freq_s=linspace(0,freq_max,t_win/2*freq_max);
 | |
|             if ibeg>1 && iend<length(c)        % eliminate indexes outside data range  
 | |
|                 [psd,integ_v,psd_int,integ_d]=integrated_psd(c(ibeg:iend,:),t_win,t_ovlp,nyqhp,dts);
 | |
|                 psd_v=psd+psd_v;
 | |
|                 psd_d=psd_int+psd_d;
 | |
|                 [frz,crsp,pwsp,coherz,nsp]=fqresp(c(ibeg:iend,:),shock_ch,t_win,t_ovlp,nyqhp,dts);
 | |
|                 frz_cut=frz+frz_cut;
 | |
|             end
 | |
|              
 | |
|         end
 | |
|         
 | |
|         psd_v=psd_v/length(t_impact);
 | |
|         psd_d=psd_d/length(t_impact);
 | |
|         frz_cut=frz_cut/length(t_impact);
 | |
|         
 | |
|      end
 | |
|     
 | |
|    
 | |
|      
 | |
|         
 | |
|     drms=max(integ_d);  % compute rms level
 | |
|     dc=hpfint(c,f_cut,dts);    % filter and integrate in time domain
 | |
|     dppc=hpdpp(dc,t_win,t_ovlp,1,dts);  % compute peak to peak level
 | |
|     
 | |
|     
 | |
|    %% transfer function, cross spectrum, power spectr. and coherence w.r.t. ch1
 | |
|       
 | |
|    eval(['c',num2str(i),'=c;'])
 | |
|    eval(['dc',num2str(i),'=dc;'])
 | |
|    eval(['dppc',num2str(i),'=dppc;'])
 | |
|    eval(['drms',num2str(i),'=drms;'])
 | |
|    eval(['psd_v',num2str(i),'=psd_v;'])  % already integrated in OROS
 | |
|    eval(['psd_d',num2str(i),'=psd_d;'])
 | |
|    eval(['integ_v',num2str(i),'=integ_v;'])
 | |
|    eval(['integ_d',num2str(i),'=integ_d;'])
 | |
|    eval(['frz',num2str(i),'=frz_cut;'])
 | |
| %   eval(['frh',num2str(i),'=frh_cut;'])
 | |
| %    eval(['frx',num2str(i),'=frx;'])
 | |
|    eval(['coherz',num2str(i),'=coherz;'])
 | |
|    eval(['time',num2str(i),'=time;'])
 | |
|    
 | |
|    clear data c dc psd psd_v psd_d time c_shut % clean up the mess
 | |
|    
 | |
| 
 | |
| end |