image thumbnail

What's New for Object-Oriented Programming in MATLAB Webinar - Code Examples

by

 

10 Apr 2008 (Updated )

Code examples used in "What's New for Object-Oriented Programming in MATLABĀ®" Webinar

handleImp.sads
% Sensor Array Data Set Handle Class
classdef sads < handle
    properties
        Data         % Sampled sensor data
        Name         % Sensor array test run name
        Spacing      % Spacing of array (m)
        SampleRate   % Sample rate (Hz)
    end
    properties (Constant)
        c=3e8;       % Speed in medium (m/s)
    end
    properties (Access=private)
        Wavelength   % Wavelength of sources (m)
    end
    properties (Dependent)
        NumSensors   % Number of sensors
        NumSamples   % Number of samples
    end
    methods
        function obj=sads(Data, Wavelength,SampleRate,Spacing,Name)
            % SADS Create sensor array data set
            % Example:
            %  sads(Data, Wavelength,SampleRate,Spacing,Name)
            obj.Data=Data;
            obj.SampleRate=SampleRate;
            obj.Spacing = Spacing;
            obj.Name=Name;
            obj.Wavelength=Wavelength;
        end
        function [mags, fflip]=magfft(obj,zeroPadTo)
            % MAGFFT Calculate the magnitude square of the FFT of the
            % sensor array sample data, zeropadding by zeroPadTo elements
            % Example:
            %  result=magfft(s,128);

            mag=zeros(obj.NumSamples,zeroPadTo); % Preallocate store of magnitudes
            f=asin((-0.5:1/(zeroPadTo-1):0.5)*obj.Wavelength/obj.Spacing)/pi; % Frequencies
            fflip=fliplr(f); % Flip frequencies

            % Take the sum over each sensor array sample
            for k=1:obj.NumSamples
                avbig=zeros(1,zeroPadTo);            % Zero pad
                avbig(1:obj.NumSensors)=obj.Data(1,:);
                response=fft(avbig)/zeroPadTo;       % FFT of normalized signal
                mag(k,:)=abs(fftshift(response)).^2; % Mag squared of FFT
            end
            mags=sum(mag);
        end
        function NumSensors=get.NumSensors(obj)
            % Get NumSensors property
            NumSensors=size(obj.Data,2);
        end
        function NumSamples=get.NumSamples(obj)
            % Get NumSamples property
            NumSamples=size(obj.Data,1);
        end
        function plot(obj)
            % PLOT Plot the sensor array sample data set
            % Example:
            %  plot(ds);
            clf
            surf(real(obj.Data)',...
                'EdgeLighting','flat','FaceLighting','none',...
                'FaceColor',[1 1 1],...
                'EdgeColor',[0 0 1]);
            view([-55.5 74]);
            xlabel('Samples')
            ylabel('Sensors')
            zlabel('Amplitude')
            xlim([1 obj.NumSamples]);
            ylim([1 obj.NumSensors]);
            set(gca,'xtick',(1:obj.NumSamples/8:obj.NumSamples));
            title(obj.Name);
        end
        function magfftplot(obj, zeroPadTo)
            % MAGFFTPLOT Plot the magnitude square of the FFT of the sensor array data
            % Example:
            %  magfftplot(s,128);
            clf
            [mags, fflip]=magfft(obj, zeroPadTo);
            semilogy(180*fflip,mags(1:zeroPadTo),'r');
            grid
            title('Averaged FFT of Sensor Data');
            xlabel('degrees');
            ylabel('amplitude');
        end
        function angles=doa(obj)
            % DOA  Estimate the direction of arrival of the sources in the
            % sensor array data set using simplistic peak finding method
            % Example:
            %  angels=doa(ds)

            zeroPadTo=256;
            [mags, fflip]=magfft(obj,zeroPadTo);
            maxtab=peakdet(mags,.1);             % Use peakdet function
            angles=sort(fflip(maxtab(:,1))*180); % Angles
        end
        function steer(obj,theta)
            % STEER Steer array electronically by angle theta (in degrees), returning a new sensor array data set
            % Example:
            %   s=steer(s,10);

            delta=obj.Spacing; thetaR=theta*(pi/180);
            Wc=2*pi*(obj.c/obj.Wavelength); % Source frequency in radians per sec
            phaseShift=exp(-j*Wc*delta*(0:obj.NumSensors-1)*sin(thetaR)/obj.c);
            obj.Data=bsxfun(@times,obj.Data,phaseShift); % Multiply by phaseshift
        end
    end
    methods (Static)
        function showarray(Targets,NumSensors,Spacing)
            % SHOWARRAY  Illustrate a sensor array with ideal sources
            % Example:
            %  showarray(Targets,NumSensors,Spacing)

            numtargs=size(Targets,1);
            orgx=(NumSensors+1)/2;

            % Plot array
            xs=(1:NumSensors)';
            ys=zeros(NumSensors,1);
            plot(xs,ys,'o');
            axis([0 NumSensors+1 0 NumSensors+1]);
            set(gca,'xtick',1:NumSensors);
            set(gca,'ytick',[]);
            hold on
            xlabel(['Sensor spacing is ',num2str(Spacing), ' m']);
            ylim([0 NumSensors/1.5])

            % Normal to array
            xs=ones(NumSensors,1).*orgx;
            ys=(0:NumSensors-1)';
            plot(xs,ys,'--k');
            title([int2str(NumSensors),' sensor array with ' int2str(numtargs) ' sources']);

            % Incident rays
            for tar=1:numtargs
                bearing=.5*pi-Targets(tar,1);
                length_line=orgx;
                ray1x=length_line*cos(bearing)+orgx;
                ray1y=length_line*sin(bearing);
                plot([ray1x orgx]', [ray1y 0]','r');
                % text(['{\it\theta}_' int2str(tar)],[ray1x ray1y])
            end

            legend('Sensors','Direction of Array','Direction of Sources');
        end
    end
end
function [maxtab, mintab]=peakdet(v, delta)
%PEAKDET Detect peaks in a vector
%        [MAXTAB, MINTAB] = PEAKDET(V, DELTA) finds the local
%        maxima and minima ("peaks") in the vector V.
%        A point is considered a maximum peak if it has the maximal
%        value, and was preceded (to the left) by a value lower by
%        DELTA. MAXTAB and MINTAB consists of two columns. Column 1
%        contains indices in V, and column 2 the found values.
% Eli Billauer, 3.4.05 (Explicitly not copyrighted).
% http://www.billauer.co.il/peakdet.html
% This function is released to the public domain; Any use is allowed.

maxtab = [];
mintab = [];

v = v(:); % Just in case this wasn't a proper vector

if (length(delta(:)))>1
    error('Input argument DELTA must be a scalar');
end

if delta <= 0
    error('Input argument DELTA must be positive');
end

mn = Inf; mx = -Inf;
mnpos = NaN; mxpos = NaN;

lookformax = 1;

for i=1:length(v)
    this = v(i);
    if this > mx, mx = this; mxpos = i; end
    if this < mn, mn = this; mnpos = i; end

    if lookformax
        if this < mx-delta
            maxtab = [maxtab ; mxpos mx];
            mn = this; mnpos = i;
            lookformax = 0;
        end
    else
        if this > mn+delta
            mintab = [mintab ; mnpos mn];
            mx = this; mxpos = i;
            lookformax = 1;
        end
    end
end
end

Contact us