%This program computes data streams in a PWM power amp with feedback clear % this is necessary to clear variables from prior runs %*************************************************************** %*************************************************************** %*** Below is the section for the user to enter and edit data. %*************************************************************** %*************************************************************** N=2^16 ;% number of time points in simulation, hence # of frequency points % N needs to be a power of 2 for the simulation to run fast. % The student version of matlab limits you to 16238 (2^14) time points. % much better accuracy will be obtained if more than 16238 time points are used % but the simulation becomes VERY slow. Vp= 10 ; % *********this is the peak (not peak-peak) output voltage of the switching stage beta=0.1 ; % *********this is the feedback factor R2/(R1+R2). % From feedback theory, the closed-loop gain will be very close to 1/beta % so the maximum input voltage you should apply should be a little less than Vp*beta. f_tri=250E3 %triangle wave frequency in Hz v_tri_mag=1 ; % magnitude of triangle wavewaveform in volts f_clock=f_tri*8 ;%clock frequency in Hz. Must be K*f_tri, where K is an integer. % fclock sets how often the simulation makes a time step % f_clock is NOT the triangle wave frequency. fclock is necessary becuase the program must % make calculations at a finite number of steps in time. So note that the in this simulation % that the pulse width modulation is now fully continuous in this simulation. If you make N larger % than 16384, you can make f_clock larger. If you make fclock larger *without* making N larger, % there will not be enough time points in one cycle of the signal, and the FFT will be noisy. loop_bandwidth=100E3% this is the feedback loop bandwidth % the feedback loop bandwidth should be as high as possible for low distortion, % but at MOST 1/5 of the triangle wave frequency. Problems arise otherwise f_unity=loop_bandwidth*v_tri_mag/Vp/beta% unity gain frequency of the 2 integrators % this relationship arises because the voltage gain of the PWM stage is (Vp/v_tri_mag) % f_zero=10E3 ; % frequency of the zero in the second integrator % *********this should be less than the feedback loop bandwidth by at least % a factor of 2, and probably by a factor of 10. f_cutoff=50E3 ; % cutoff frequency of the output filter % usually set to somewhat higher than the maximum input signal frequency n_signal_cycles=16*8; % the signal frequency is given by f_signal=n_signal_cycles/tfinal % specifying the signal frequency in this fashion forces an interger # of % signal periods in the fft, which is necessary to obtain a clean fft. % So, if you change N or f_tri, you will need to change n_signal_cycles in order % to get the signal frequency you want. v_in = 0.8 ;% the input signal voltage amplitude, peak, not RMS gm1=2E-3 ; %transconductance of the first Gm element gm2=2E-3 ; %transconductance of the second Gm element c1=gm1*0.159/f_unity % first integration capacitor Rzero=1/gm2 ; % the resistor which introduces a zero in the second integrator % this forces the 2nd integrator gain to unity at frequencies *above* the zero c2=0.159/f_zero/Rzero % second integration capacitor Idac1_val=beta*gm1*Vp ; % % If the feedback is in the form of a voltage divider of % attenuation beta, feeding back to the *input* of the first gm stage, % then Idac1_val should be given by Idac1_val=beta*gm1*Vp % This is the configuration most probably used by 137b students %*************************************************************** %*************************************************************** %*************************************************************** %*************************************************************** %The program now sets up time and frequency increments for the simulation. f_nyquist=2*f_clock ; fhigh=2*f_nyquist ; tfinal=1/fhigh*(N-1)^2/N ;% f_signal=n_signal_cycles/tfinal % timestep=tfinal/(N-1); % timestep f=linspace(0,fhigh,N); delta_f=1/tfinal % frequency step t=linspace(0,tfinal,N); N_clock=f_clock*tfinal ; % number of data bits in simulation t_clock=1/f_clock ; % the clock period. %************************************************************************* %************************************************************************* % the values of these variables are first set to very small values % before running the finite-time simulations. %Seems to be necessary to avoid a progam bug. % The real values are subsequently calculated in the next block down v1=1E-9*sin(2*pi*f_clock*t); % v1 is the voltage on the first integration capacitor C1 v2=1E-9*sin(2*pi*f_clock*t); % v2 is the voltage on the second integration capacitor C2 v3=1E-9*sin(2*pi*f_clock*t); % v3 is the voltage at the output of the second integrator v_tri=1E-9*sin(2*pi*f_clock*t); % v_clock is the triange wave generator v_data=1E-9*sin(2*pi*f_clock*t); v_data2=1E-9*sin(2*pi*f_clock*t); %v_data(n_count)=1 ; % v_data is the output of the latched comparator Idac1=1E-9*sin(2*pi*f_clock*t); % the output of the first DAC Idac2=1E-9*sin(2*pi*f_clock*t); % the output of the second DAC %************************************************************************* %************************************************************************* % now we simulate the feedback loop in finite time steps % this is the core of the simulation % First apply an input signal: a sinewave input_signal=v_in*sin(2*pi*f_signal*t); % v_tri=v_tri_mag*( sin(2*pi*f_tri*t)+sin(4*pi*f_tri*t) ); % time=0 ; n_count=1; % intializing the count variable time_int=0; % this is accumulated time within a clock cycle...see below % we now have a if/else/for loop which increments time in small steps for n_count=2:N ; squarewave=1 ; n_countless=n_count-1; time=n_count*timestep; % increment time by one timestep input_signal(n_count)=v_in*sin(2*pi*f_signal*time);% the input signal v_tri(n_count)=cos(2*pi*f_tri*time); % the triangle wave: fundamental v_tri(n_count)=v_tri(n_count)+(1/9)*cos(3*2*pi*f_tri*time); % the triangle wave: adding the 3rd harmonic v_tri(n_count)=v_tri(n_count)+(1/25)*cos(5*2*pi*f_tri*time); % the triangle wave: adding the 5th harmonic v_tri(n_count)=v_tri(n_count)+(1/49)*cos(7*2*pi*f_tri*time); % the triangle wave: adding the 7th harmonic v_tri(n_count)=v_tri(n_count)+(1/81)*cos(9*2*pi*f_tri*time); % the triangle wave: adding the 9th harmonic v_tri(n_count)=v_tri(n_count)+(1/121)*cos(11*2*pi*f_tri*time); % the triangle wave: adding the 11th harmonic v_tri(n_count)=v_tri(n_count)+(1/169)*cos(13*2*pi*f_tri*time); % the triangle wave: adding the 13th harmonic v_tri(n_count)=v_tri(n_count)+(1/225)*cos(15*2*pi*f_tri*time); % the triangle wave: adding the 15th harmonic v_tri(n_count)=v_tri(n_count)+(1/289)*cos(17*2*pi*f_tri*time); % the triangle wave: adding the 17th harmonic v_tri(n_count)=v_tri_mag*v_tri(n_count)/1.2; % the triangle wave: scaling for the correct amplitude if v_data(n_countless)> 0 ; % setting feedback signal according to the comparator output Idac1(n_count)=Idac1_val; else Idac1(n_count)=-1*Idac1_val; end % summing currents at C1, and then doing I = c dV/dt to find the voltage on c1 I_nodeone(n_count)= gm1*input_signal(n_count)-Idac1(n_count); n_count_less=n_count-1; v1(n_count)=v1(n_count_less)+I_nodeone(n_count)*timestep/c1; % finding the current on C2 and then doing I = c dV/dt to find the voltage on c2 I_nodetwo(n_count)= gm2*v1(n_count); v2(n_count)=v2(n_count_less)+I_nodetwo(n_count)*timestep/c2; % adding the voltage drop on the resistor Rzero v3(n_count)=v2(n_count)+Rzero*gm2*v1(n_count) ; time_int=time_int+timestep; % time accumulated within the clock period if time_int < t_clock ;% if the accumulated time within a bit is less than a clock period v_data(n_count)=v_data(n_count_less) ;% then the data keeps its same value v_data2(n_count)=squarewave*v_data(n_count); else % otherwise the data changes according to the new comparator input v_data(n_count)=Vp*sign(v3(n_count)-v_tri(n_countless));%this line for 2 integrators %v_data(n_count)=Vp*sign(v1(n_count)-v_tri(n_countless));% this line for 1 integrator time_int = time_int -t_clock ; % and the accumlated time is re-incremented (new clock!) end end %************************************************************************* %************************************************************************* % now we FFT the data stream, and then filter it so that we can see the % resulting low-pass-filtered output signal % % first we must shuffle the frequencies in the modulation filter % to get them in the strange order the FFT requires alpha=0.51 ; n_count=1; for n_count = 1:N if n_count > N/2 f_c=1/tfinal*(N-1)/N*(n_count-N); else f_c=1/tfinal*(N-1)/N*(n_count-1); end H_filter(n_count)=raised_cos_baseband_filter2(f_c,f_cutoff,alpha); % root-raised-cos filter % the filter I have used for convenience is a mathematically idealized raised-cosine filter. % your real circuit has the low pass filter determined by the values of L and C you use % on the output. end % %************************************************************************* %now we can filter the modulated signal v_data_f=fft(v_data) ;% first FFT the data stream input_signal_f=fft(input_signal); % also handy to have the fft of the input signal % % n_count=1; for n_count=1:N % then multiply it with the filter transfer function. filt_data_f(n_count)=v_data_f(n_count)*H_filter(n_count) ; end filt_data_t=ifft(filt_data_f); % and then IFFT to get the time waveform % % in_spec = fft(input_signal.*hann(N)); out_spec = fft(v_data.*hann(N)); % % %************************************************************************* %************************************************************************* %************************************************************************* %************************************************************************* %************************************************************************* % %Plotting calls are below % close all %*************************************************************** % **************** first plot window *************************** scrsz = get(0,'ScreenSize'); %figure('Position',[1 scrsz(4)/2 scrsz(3)/2 scrsz(4)/2]) % left bottom width height fg1=figure('Position',[0 0 scrsz(3)/1.1 scrsz(4)/1.1]); % % this plots the first integrator output voltage ndatasplot=30; % enter here the number of signal cyles you want plotting t_sigplot_stop=ndatasplot*t_clock; plot(t,v1) xlabel('time, sec') %axis([0 t_sigplot_stop -2.5 2.5]) %axis([0 tfinal -2 2]) title('first integrator output voltage: MAKE SURE IT IS LESS THAN THE CIRCUIT Voltage CLIPPING LIMITS') % %*************************************************************** % **************** 2nd plot window *************************** scrsz = get(0,'ScreenSize'); %figure('Position',[1 scrsz(4)/2 scrsz(3)/2 scrsz(4)/2]) % left bottom width height fg2=figure('Position',[0 0 scrsz(3)/1.1 scrsz(4)/1.1]); % this plots the output voltage of the 2nd integrator ndatasplot=30; % enter here the number of signal cyles you want plotting t_sigplot_stop=ndatasplot*t_clock; plot(t,v3) xlabel('time, sec') %axis([0 t_sigplot_stop -2.5 2.5]) %axis([0 tfinal -2 2]) title('second integrator output voltage: MAKE SURE IT IS LESS THAN THE CIRCUIT Voltage CLIPPING LIMITS') % %*************************************************************** % **************** 3rdplot window *************************** scrsz = get(0,'ScreenSize'); %figure('Position',[1 scrsz(4)/2 scrsz(3)/2 scrsz(4)/2]) % left bottom width height fg3=figure('Position',[0 0 scrsz(3)/1.1 scrsz(4)/1.1]); % % this plots the ouput voltage unfiltered and filtered ndatasplot=600; % number of data cyles you want plotting t_sigplot_stop=ndatasplot*t_clock; plot(t,v_data,t,filt_data_t) xlabel('time, sec') axis([0 t_sigplot_stop -1.1*Vp 1.1*Vp]) %axis([0 tfinal -2 2]) title('ouput voltage, unfiltered digtal swiching pattern and filtered analog output voltage') % %*************************************************************** % **************** 4th plot window *************************** scrsz = get(0,'ScreenSize'); %figure('Position',[1 scrsz(4)/2 scrsz(3)/2 scrsz(4)/2]) % left bottom width height fg4=figure('Position',[0 0 scrsz(3)/1.1 scrsz(4)/1.1]); % % this plots input signal voltage ndatasplot=30; % enter here the number of signal cyles you want plotting t_sigplot_stop=ndatasplot*t_clock; plot(t,input_signal) xlabel('time, sec') %axis([0 t_sigplot_stop -2.5 2.5]) %axis([0 tfinal -2 2]) title('input signal voltage ') % % %*************************************************************** %*************************************************************** % **************** 5th plot window *************************** fg5=figure('Position',[0 0 scrsz(3)/1.1 scrsz(4)/1.1]); %******************************************** % this plots the spectrum of the data %semilogx( f ,20*log10(2*abs(v_data_f)/N+1e-90),f, ... % 20*log10(2*abs(input_signal_f)/N+1e-90)); semilogx( f ,20*log10(2*abs(in_spec)/N+1e-90),f, ... 20*log10(2*abs(out_spec)/N+1e-90)); %axis([ 0 5*f_clock -50 -10 ]) %axis([ 0 log10(fhigh/2) -100 -10 ]) axis([ 0 fhigh/2 -100 20*log10(abs(Vp)) ]) %set('XScale', 'log) xlabel('Frequency, Hz') ylabel('Power, dB_V , with delta-f Hz FFT frequency resolution ') % you need to work out the FFT frequency resolution !!!!! title('input and output signals in the frequency domain') % %*************************************************************** %*************************************************************** % **************** 6th plot window *************************** scrsz = get(0,'ScreenSize'); %figure('Position',[1 scrsz(4)/2 scrsz(3)/2 scrsz(4)/2]) % left bottom width height fg6=figure('Position',[0 0 scrsz(3)/1.1 scrsz(4)/1.1]); % % ndatasplot=10; % enter here the number of signal cyles you want plotting t_sigplot_stop=ndatasplot*t_clock; plot(t,I_nodeone) xlabel('time, sec') %axis([0 t_sigplot_stop -2.5 2.5]) %axis([0 tfinal -2 2]) title('first gm stage output current: MAKE SURE IT IS LESS THAN THE CIRCUIT CURRENT CLIPPING LIMITS') %*************************************************************** % **************** 7th plot window *************************** scrsz = get(0,'ScreenSize'); %figure('Position',[1 scrsz(4)/2 scrsz(3)/2 scrsz(4)/2]) % left bottom width height fg7=figure('Position',[0 0 scrsz(3)/1.1 scrsz(4)/1.1]); % % this plots second gm stage output current ndatasplot=10; % enter here the number of signal cyles you want plotting t_sigplot_stop=ndatasplot*t_clock; plot(t,I_nodetwo) xlabel('time, sec') %axis([0 t_sigplot_stop -2.5 2.5]) %axis([0 tfinal -2 2]) title('second gm stage output current: MAKE SURE IT IS LESS THAN THE CIRCUIT CURRENT CLIPPING LIMITS') % %*************************************************************** % **************** 8th plot window *************************** scrsz = get(0,'ScreenSize'); %figure('Position',[1 scrsz(4)/2 scrsz(3)/2 scrsz(4)/2]) % left bottom width height fg8=figure('Position',[0 0 scrsz(3)/1.1 scrsz(4)/1.1]); % % this plots second gm stage output current ndatasplot=150; % enter here the number of signal cyles you want plotting t_sigplot_stop=ndatasplot*t_clock; plot(t,v_tri,t,v1,t,v_data/Vp ) xlabel('time, sec') axis([0 t_sigplot_stop -2.5 2.5]) %axis([0 tfinal -2 2]) title('triangle wave (blue), comparator input (green), normalized comparator output (red)') %