fseriesexpogui.m

非常好的数字处理教程

function varargout = fseriesexpogui(varargin)
% This code is by Jim Squire <[e-mail]>
%	It needs to be recoded to SSUM format.
% 
% FSERIESEXPOGUI M-file for fseriesexpogui.fig
%      FSERIESEXPOGUI, by itself, creates a new FSERIESEXPOGUI or raises the existing
%      singleton*.
%
%      H = FSERIESEXPOGUI returns the handle to a new FSERIESEXPOGUI or the handle to
%      the existing singleton*.
%
%      FSERIESEXPOGUI('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in FSERIESEXPOGUI.M with the given input arguments.
%
%      FSERIESEXPOGUI('Property','Value',...) creates a new FSERIESEXPOGUI or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before fseriesexpogui_OpeningFunction gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to fseriesexpogui_OpeningFcn via varargin.
%
%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one
%      instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help fseriesexpogui

% Last Modified by GUIDE v2.5 31-Mar-2005 18:06:31

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
    'gui_Singleton',  gui_Singleton, ...
    'gui_OpeningFcn', @fseriesexpogui_OpeningFcn, ...
    'gui_OutputFcn',  @fseriesexpogui_OutputFcn, ...
    'gui_LayoutFcn',  @fseriesexpogui_LayoutFcn, ...
    'gui_Callback',   []);
if nargin & isstr(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT


% --- Executes just before fseriesexpogui is made visible.
function fseriesexpogui_OpeningFcn(hObject, eventdata, handles, varargin)
handles.output = hObject;
guidata(hObject, handles);
t=load('fseriespict.mat');
% Fig1
axes(handles.uiaxFig1);
image(t.strGraphic.Fig1);
set(handles.uiaxFig1,'Visible','Off')
colormap([0 0 0;236 233 216]/255)
% Fig2
axes(handles.uiaxFig2)
image(t.strGraphic.Fig2)
set(handles.uiaxFig2,'Visible','Off')
% Fig3
axes(handles.uiaxFig3)
image(t.strGraphic.Fig3)
set(handles.uiaxFig3,'Visible','Off')
% Calculate and fill the plots
CalculateAndUpdate(handles)

function varargout = fseriesexpogui_OutputFcn(hObject, eventdata, handles)
varargout{1} = handles.output;


%%%%%%%%%%%%%%%%%%%%%%%% Input function: wave type %%%%%%%%%%%%%%%%%%%%%%%%%%%
function uirbWaveTypeSin_Callback(hObject, eventdata, handles)
set(handles.uirbWaveTypeSin, 'Value', 1);
set([handles.uirbWaveTypeSquare handles.uirbWaveTypeTriangle], 'Value', 0);
CalculateAndUpdate(handles)

function uirbWaveTypeSquare_Callback(hObject, eventdata, handles)
set(handles.uirbWaveTypeSquare, 'Value', 1);
set([handles.uirbWaveTypeSin handles.uirbWaveTypeTriangle], 'Value', 0);
CalculateAndUpdate(handles)

function uirbWaveTypeTriangle_Callback(hObject, eventdata, handles)
set(handles.uirbWaveTypeTriangle, 'Value', 1);
set([handles.uirbWaveTypeSin handles.uirbWaveTypeSquare], 'Value', 0);
CalculateAndUpdate(handles)

%%%%%%%%%%%%%%%%%%%%%%% Slider data %%%%%%%%%%%%%%%%%%%%%%%%%
function uislFrequency_Callback(hObject, eventdata, handles)
CalculateAndUpdate(handles)

function uislAmplitude_Callback(hObject, eventdata, handles)
CalculateAndUpdate(handles)

function uislSymmetry_Callback(hObject, eventdata, handles)
CalculateAndUpdate(handles)

function uislYOffset_Callback(hObject, eventdata, handles)
CalculateAndUpdate(handles)

function uislTDelay_Callback(hObject, eventdata, handles)
CalculateAndUpdate(handles)

%%%%%%%%%%%%%%%%%%%%%%% Fourier series form %%%%%%%%%%%%%%%%%%%%%%%%%
function uirbFig1_Callback(hObject, eventdata, handles)
set(handles.uirbFig1, 'Value', 1)
set([handles.uirbFig2 handles.uirbFig3], 'Value', 0);
set(handles.uirbPlotComponent1,'String','ai')
set(handles.uirbPlotComponent2,'String','bi')
set([handles.uirbPlotComponent3 handles.uirbPlotComponent4],'Enable','off')
if get(handles.uirbPlotComponent3,'Value') ||get(handles.uirbPlotComponent4,'Value') 
    set([handles.uirbPlotComponent3 handles.uirbPlotComponent4],'Value',0)
    set(handles.uirbPlotComponent1,'Value',1)
end
CalculateAndUpdate(handles)

function uirbFig2_Callback(hObject, eventdata, handles)
set(handles.uirbFig2, 'Value', 1);
set([handles.uirbFig1 handles.uirbFig3], 'Value', 0);
set(handles.uirbPlotComponent1,'String','Ai')
set(handles.uirbPlotComponent2,'String','Thetai')
set([handles.uirbPlotComponent3 handles.uirbPlotComponent4],'Enable','off')
if get(handles.uirbPlotComponent3,'Value') ||get(handles.uirbPlotComponent4,'Value') 
    set([handles.uirbPlotComponent3 handles.uirbPlotComponent4],'Value',0)
    set(handles.uirbPlotComponent1,'Value',1)
end
CalculateAndUpdate(handles)

function uirbFig3_Callback(hObject, eventdata, handles)
set(handles.uirbFig3, 'Value', 1);
set([handles.uirbFig1 handles.uirbFig2], 'Value', 0);
set(handles.uirbPlotComponent1,'String','real{ci}')
set(handles.uirbPlotComponent2,'String','imag{ci}')
set([handles.uirbPlotComponent3 handles.uirbPlotComponent4],'Enable','on')
CalculateAndUpdate(handles)


%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Fourier Series: NumTerms %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function uislNumTerms_Callback(hObject, eventdata, handles)
switch round(get(handles.uislNumTerms,'Value'))
    case 0, N = 0;
    case 1, N = 1;
    case 2, N = 2;
    case 3, N = 3;
    case 4, N = 4;
    case 5, N = 5;
    case 6, N = 6;
    case 7, N = 7;
    case 8, N = 8;
    case 9, N = 9;
    case 10, N = 10;
    case 11, N = 20;
    case 12, N = 30;
    case 13, N = 40;
    case 14, N = 50;
    case 15, N = 100;
end
set(handles.uitxNumTerms,'String',sprintf('N terms = %g',N))
CalculateAndUpdate(handles)

%%%%%%%%%%%%%%%%%%%%%%%%%% Fourier Series: plot components %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function uirbPlotComponent1_Callback(hObject, eventdata, handles)
set(handles.uirbPlotComponent1, 'Value', 1);
set([handles.uirbPlotComponent2 handles.uirbPlotComponent3 handles.uirbPlotComponent4], 'Value', 0);
CalculateAndUpdate(handles)

function uirbPlotComponent2_Callback(hObject, eventdata, handles)
set(handles.uirbPlotComponent2, 'Value', 1);
set([handles.uirbPlotComponent1 handles.uirbPlotComponent3 handles.uirbPlotComponent4], 'Value', 0);
CalculateAndUpdate(handles)

function uirbPlotComponent3_Callback(hObject, eventdata, handles)
set(handles.uirbPlotComponent3, 'Value', 1);
set([handles.uirbPlotComponent1 handles.uirbPlotComponent2 handles.uirbPlotComponent4], 'Value', 0);
CalculateAndUpdate(handles)

function uirbPlotComponent4_Callback(hObject, eventdata, handles)
set(handles.uirbPlotComponent4, 'Value', 1);
set([handles.uirbPlotComponent1 handles.uirbPlotComponent2 handles.uirbPlotComponent3], 'Value', 0);
CalculateAndUpdate(handles)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% CalculateAndUpdate %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function CalculateAndUpdate(handles)
strGui = GetGUIValues(handles);  % get all the scaled date from the GUI

% calculate strData.y0Time, .xTime, y1Time, xFreq, yFreq
strData = CaculateSeries(strGui);

% Display time
axes(handles.uiaxTime)
plot(strData.xTime,strData.y0Time,'r-',strData.xTime,strData.y1Time,'b-')
xlabel('Time (s)')
legend('Orignal','Reconstructed')
v=axis;
if v(4)-v(3)<1e-4
    v(4)=v(4)+0.5;
    v(3)=v(3)-0.5;
end
axis(v)
% Display freq
axes(handles.uiaxFreq)
stem('v6',strData.xFreq,strData.yFreq)
xlabel('Freq (Hz)')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% GetGUIValues %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function strGui = GetGUIValues(handles)
% Wave type
t = get(handles.uirbWaveTypeSin,'Value') + 2*get(handles.uirbWaveTypeSquare,'Value') + 4*get(handles.uirbWaveTypeTriangle,'Value');
switch t
    case 1, strGui.waveType = 'Sin';
    case 2, strGui.waveType = 'Square';
    case 4, strGui.waveType = 'Triangle';
    otherwise, error(sprintf('Unknown wave type: ',t));
end
% Wave parameters
strGui.w = ScaleIt(get(handles.uislFrequency,'Value'),0,4,8,0,2*pi,6*pi);                  % frequency in rads/s
if ~strGui.w, strGui.w=1e-6; end
T = 2*pi/strGui.w;
strGui.A = ScaleIt(get(handles.uislAmplitude,'Value'),0,4,8,0,1,4);                        % amplitude
strGui.sym = ScaleIt(round(get(handles.uislSymmetry,'Value')),0,4,8,1e-6,0.5-1e-6,1-1e-6); % symmetry
strGui.yOff = ScaleIt(get(handles.uislYOffset,'Value'),0,4,8,-2,0,2);                      % yOffset
strGui.tDelay = ScaleIt(get(handles.uislTDelay,'Value'),0,4,8,0,T/2,T);                    % tDelay
% Fourier series form
t = get(handles.uirbFig1,'Value') + 2*get(handles.uirbFig2,'Value') + 4*get(handles.uirbFig3,'Value');
switch t
    case 1, strGui.fourierForm = 'Quadrature';
    case 2, strGui.fourierForm = 'Polar';
    case 4, strGui.fourierForm = 'Exponential';
    otherwise, error(sprintf('Unknown Fourier Form',t));
end
% Num Terms
sNTerms = get(handles.uitxNumTerms,'String');
sNTerms(1:10) = [];
strGui.N = str2num(sNTerms);
% Plot component ai, bi, etc.
t = get(handles.uirbPlotComponent1,'Value') + 2*get(handles.uirbPlotComponent2,'Value') + 4*get(handles.uirbPlotComponent3,'Value') + 8*get(handles.uirbPlotComponent4,'Value');
switch t
    case 1, strGui.plotComponent = 1;
    case 2, strGui.plotComponent = 2;
    case 4, strGui.plotComponent = 3;
    case 8, strGui.plotComponent = 4;
    otherwise, error(sprintf('Unknown component type %g',t))
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% CaculateSeries %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function strData = CaculateSeries(strGui)

% get GUI parameters
waveType = strGui.waveType;
w = strGui.w;
T = 2*pi/w;
A = strGui.A;
sym = strGui.sym;
yOff = strGui.yOff;
tDelay = strGui.tDelay;
fourierForm = strGui.fourierForm;
N = strGui.N;
plotComponent = strGui.plotComponent;

%%%%%%%%%%%%%%%%%%%%%%%%%%% Calculate the time-domain vectors %%%%%%%%%%%%%%%%%%%%%
% create output vectors
t = linspace(0,2,400);
y0 = zeros(size(t));
y1 = y0;
% First create a triangle wave, regardless of what was requested
T1 = T*sym;  % Triangle changes direction at T1, which may be from 0 to T
m1 = -2/T1; b1 = 1; % line equations for first part
m2 = 2/(T-T1); b2 = 1-m2*T;
tt = mod(t-tDelay,T) + tDelay;
for i=1:length(t)
    if tDelay<=tt(i) && tt(i)<=T1+tDelay
        triangle(i) = m1*(tt(i)-tDelay) + b1;
    elseif T1+tDelay < tt(i) && tt(i) <= T+tDelay
        triangle(i) = m2*(tt(i)-tDelay) + b2;
    else
        error('Impossible time bracket')
    end
end
% Next change it into the requested type
switch strGui.waveType
    case 'Sin'
        for i=1:length(y0)
            theta = (1-triangle)*pi/2;
            y0 = A*cos(theta) +yOff;
        end
    case 'Square'
        if sym<=1e-6,
            y0 = ones(size(triangle));
        elseif sym>=1-1e-6
            y0 = -ones(size(triangle));
        else
            symOffset = (sym-1/2)*2;
            y0 = A*sign(triangle-symOffset) + yOff;
        end
        Ta = (m1*tDelay-b1+2*sym-1)/m1;
        Tb = (m2*tDelay-b2+2*sym-1)/m2;
        if w<=1e-6
            Ta=0; Tb=0;
        end
    case 'Triangle'
        y0 = A*triangle + yOff;
    otherwise
        error('Unknown waveType')
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Calculate the c's %%%%%%%%%%%%%%%%%%%%%%%%%%%%
% create vectors
f = zeros(1,N+1);
c = f;
% determine c
switch strGui.waveType
    case 'Sin'
        tt1=m1*m2*pi*T*yOff;
        tt2=m2*cos(b1*pi/2)-m1*cos(pi/2*(b2+m2*T));
        tt3=-m2*cos(pi/2*(b1+m1*T1))+m1*cos(pi/2*(b2+m2*T1));
        tt4=m1*m2*pi*T;
        c(1)=(tt1+2*A*(tt2+tt3))/tt4;
        if w<=1e-5
            c(1)=y0(1);
        else
            for n=1:N
                tt1=-1/(n*T*w);
                tt2=pi/2;
                tt3=m1^2*pi^2-4*n^2*w^2;
                tt4=m2^2*pi^2-4*n^2*w^2;
                tt5=tt2*(b1+m1*T1);
                tt6=tt2*(b2+m2*T);
                tt7=tt2*(b2+m2*T1);
                tt8=j*m1*pi*cos(b1*tt2)-2*n*w*sin(b1*tt2);
                tt9=-j*m1*pi*cos(tt5)+exp(j*n*T1*w)*tt8+2*n*w*sin(tt5);
                tt10=-j*m2*pi*cos(tt6)+2*n*w*sin(tt6);
                tt11=j*m2*pi*tt3*cos(tt7)+tt4*tt9-2*n*w*tt3*sin(tt7);
                tt12=exp(-j*n*w*(T+tDelay))*yOff*(exp(j*n*T*w)-1);
                tt13=exp(j*n*w*(T1+2*tDelay))*tt3*tt10+exp(j*n*w*(T+2*tDelay))*tt11;
                tt14=2*A*n*w*tt13*exp(-j*n*w*(T+T1+3*tDelay));
                c(n+1)=tt1*j*(tt14/(tt3*tt4));
            end
        end
    case 'Triangle'
        if w<=1e-5
            c(1)=y0(1);
        else
            c(1)=A*(m2*T^2+2*b2*(T-T1)+2*b1*T1+(m1-m2)*T1^2)/(2*T)+yOff;
            for n=1:N
                tt1=j*n*w;
                tt2=exp(-tt1*(T+T1+tDelay))/(n^2*T*w^2);
                tt3=-exp(tt1*(T+T1))*(m1+tt1*b1);
                tt4=exp(tt1*T1)*(m2+tt1*(b2+m2*T));
                tt5=exp(tt1*T)*(m1-m2+tt1*(b1-b2+m1*T1-m2*T1));
                tt6=-tt1*exp(tt1*T1)*(exp(tt1*T)-1)*yOff;
                c(n+1)=tt2*(A*(tt3+tt4+tt5)+tt6);
            end
        end
    case 'Square'
        c(1)=A*(T+2*(Ta-Tb))/T+yOff; % this is really c0
        for n=1:N
            tt1=j*n*w;
            tt2=j*exp(-tt1*(T+2*Ta+Tb))/(n*T*w);
            tt3=-2*A*exp(tt1*(T+2*Ta));
            tt4=A*exp(tt1*(Ta+Tb));
            tt5=(A-yOff)*exp(tt1*(T+Ta+Tb));
            tt6=yOff*exp(tt1*(Ta+Tb));
            c(n+1)=tt2*(tt3+tt4+tt5+tt6);
        end
    otherwise
        error('unknown waveType')
end
% if elements approximately zero, make exactly zero
for n=1:N+1
    if abs(real(c(n)))<=1e-5
        c(n)=j*imag(c(n));
    end
    if abs(imag(c(n)))<=1e-5
        c(n)=real(c(n));
    end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Synthesize f(t) %%%%%%%%%%%%%%%%%%%%%%%%%%%%
y1=y1+c(1);
for n=1:N
    y1=y1 + c(n+1)*exp(j*n*w*t) + conj(c(n+1))*exp(-j*n*w*t);
end
switch fourierForm
    case 'Quadrature'
        switch plotComponent
            case 1 %an
                yFreq=2*real(c);
            case 2 %bn
                yFreq=2*imag(c);
                yFreq(1)=0;
            otherwise
                error('bad plotComponent combination')
        end
    case 'Polar'
        switch plotComponent
            case 1 %An
                yFreq=2*abs(c);
            case 2 %theta
                yFreq=angle(c)*180/pi;
            otherwise
                error('bad plotComponent combination')
        end
    case 'Exponential'
        switch plotComponent
            case 1 %real
                yFreq=real(c);
            case 2 %imag
                yFreq=imag(c);
            case 3 %mag
                yFreq=abs(c);
            case 4 %angle
                yFreq=angle(c)*180/pi;
            otherwise
                error('bad plotComponent combination')
        end
end
if w<=1e-6
    xFreq=zeros(size(c));
else
    xFreq= linspace(0,N*w/(2*pi),N+1);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%% Fill up output structure %%%%%%%%%%%%%%%%%%%%%%
strData.xTime = t;
strData.y0Time = y0;
strData.y1Time = y1;
strData.xFreq = xFreq;
strData.yFreq = yFreq;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ScaleIt %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function out=ScaleIt(sliderIn,sliderMin,sliderMean,sliderMax,outMin,outMean,outMax)
% ScaleIt provides a piecewise scaling of value SliderIn
% if SliderIn is between sliderMin and SliderMean it is scaled between outMin and outMean
% and similarly betwen Mean and Max
if sliderIn<sliderMin || sliderIn > sliderMax
    error('slider in outside min/max boundaries in ScaleIt')
end
if sliderIn < sliderMean
    tm = (outMean-outMin)/(sliderMean-sliderMin);
else
    tm = (outMax-outMean)/(sliderMax-sliderMean);
end
tb = outMean - (sliderMean * tm);
out = sliderIn * tm + tb;



% --- Creates and returns a handle to the GUI figure. 
function h1 = fseriesexpogui_LayoutFcn(policy)
% policy - create a new figure or use a singleton. 'new' or 'reuse'.

persistent hsingleton;
if strcmpi(policy, 'reuse') & ishandle(hsingleton)
    h1 = hsingleton;
    return;
end

appdata = [];
appdata.GUIDEOptions = struct(...
    'active_h', [], ...
    'taginfo', struct(...
    'figure', 2, ...
    'text', 34, ...