Explanation of the matlab code

Hi I foun this 'extract' of a matlab code on the internet. Can someone explain to me line by line whats happening? Im REALLY in need of help.

*

% ------------- % This is code to make the edge detecting filter % ----%

function filter=gaussfilt(N)

% calculate alpha so the filter fills all N points

alpha=N;

first=-(1-N/2)*exp(-(1-N/2)^2/alpha);

count=0;

while first<.1*(-(1530/4000*N-N/2)*exp(-(1530/4000*N-N/2)^2/alpha))

      count=count+1;

      alpha=N*500*count;

      first=-(1-N/2)*exp(-(1-N/2)^2/alpha);

end

for n=1:N

      filter(n)=-(n-N/2)*exp(-(n-N/2)^2/alpha);   % d/dt of a gaussian

end

filter=filter/sum(abs(filter));     % normalization

return

[John Doe]

http://stackoverflow.com/a/17134786/2338750

[Jan]

@Jolini: Please take the time to ask specific questions. I'm sure it is a waste of time to explain the line "alpha = N;".

[Honglei Chen]

Wherever you get the code, it is a bad idea that you use filter as the variable name since this is also a built-in function name.

[Jolini]

Yes, i realized that. I changed the variable name. Thank you.

[Muthu Annamalai]

This not mandatory, MATLAB can understand the scoping of variables, and the shadowing of functions. You are free to use symbol names of your choice - but - metasyntactic meanings are useful to other users of your code, and maintainers. After all a computer can understand stuff in binary!

[Muthu Annamalai]

Also I would recommend you use code from MATLAB File Exchange, and share it back. Its a great community, and you may like it too! http://www.mathworks.com/matlabcentral/fileexchange/

[Md Ibrahim]

vmax=8;
vmin=-vmax;
del=(vmax-vmin)/L;
part=vmin:del:vmax;    % level are between vmin and vmax with difference of del
code=vmin-(del/2):del:vmax+(del/2);        % Contain Quantized values 
[ind,q]=quantiz(s,part,code);                     % Quantization process
% ind contain index number and q contain quantized  values
l1=length(ind);
l2=length(q);

for i=1:l1
    if(ind(i)~=0)       % To make index as binary decimal so started from 0 to N
        ind(i)=ind(i)-1;
    end 
    i=i+1;
end   
for i=1:l2
    if(q(i)==vmin-(del/2))       % To make quantize value in between the levels
        q(i)=vmin+(del/2);
    end
end    
subplot(3,1,3);
stem(q);grid on;                                 % Display the Quantize values
title('Quantized Signal');
ylabel('Amplitude--->');
xlabel('Time--->');

[Murthy]

% Load the voice signal [voice, fs] = audioread('stomp.mp3');

% Add sinusoidal interference at 1500Hz t = (0:length(voice)-1)/fs; sin_interference = 0.1*sin(2*pi*1500*t); corrupted_voice = voice , sin_interference;

% Design FIR null filter n = 100; % Filter length f = [0 0.7 0.75 1]; % Normalized frequency bands a = [1 1 0 0]; % Amplitude at each band b = firpm(n, f, a);

% Filter the corrupted signal to recover the original voice signal recovered_voice = filter(b, 1, corrupted_voice);

% Play the original, corrupted, and recovered voice signals sound(voice, fs); pause(1); sound(corrupted_voice, fs); pause(10); sound(recovered_voice, fs); pause(50) % Plot the original, corrupted, and recovered voice signals figure; subplot(3,1,1); plot(voice); title('Original Voice Signal'); subplot(3,1,2); plot(corrupted_voice); title('Corrupted Voice Signal'); subplot(3,1,3); plot(recovered_voice); title('Recovered Voice Signal');

% Plot the frequency response of the null filter figure; freqz(b, 1); title('Frequency Response of Null Filter');

Answer 1

[Muthu Annamalai]

Your code calculates a filter kernel. I haven't run the code, but, it seems to me a kind of Gaussian filter, with some normalizations.

It could be rewritten like,

    function filter=gaussfilt(N)

    % calculate alpha so the filter fills all N points

    alpha=N;

    first=-(1-N/2)*exp(-(1-N/2)^2/alpha);

    count=0;
    % find the standard deviation for the Gaussian
    while first<.1*(-(1530/4000*N-N/2)*exp(-(1530/4000*N-N/2)^2/alpha))

          count=count+1;

          alpha=N*500*count;

          first=-(1-N/2)*exp(-(1-N/2)^2/alpha);

    end

    % calculate the filter kernel
    n=1:N
    filter =-(n-N/2).*exp(-(n-N/2).^2/alpha); 
    %normalize kernel & return

HTH

[Jolini]

Hi. Thank you for the reply. Yes, its a Gaussian filter. I anyhow made some changes and I've finally understood it. I'm having trouble in finding peaks (maximum) of an enveloped signal. I've posted the question in MATLAB.. Would you be able to go through it and possibly help me out??

This is the link.

http://www.mathworks.us/matlabcentral/answers/79586-finding-peaks-of-a-signal

Thanx in advance

[Usama]

H1_y= (I./(2*pi.*((X-r*cos((n-1)*2*pi/N)).^2+(Y-r*sin((n-1)*2*pi/N)).^2))).*(X-r*cos((n-1)*2*pi/N)).*...
    (cos((n-1)*2*pi/N));

Answer 2

[varshini rebala]

Can any one explain me what is happing in this code ?

% Selected Mapping (SLM)is one of the techniques used for % peak-to-average power ratio (PAPR) reduction in OFDM systems. % In this .m file it is simulated for a QPSK modulated, 64-subband % OFDM symbols.

clc clear

load ofdm_100000 % this is a .mat file containing 100000 QPSK modulated, % 64-element OFDM symbols. It is available with the % partial_transmit_sequence.m file, previously submitted % by the auther, for free download at the 'file exchange' % web page.

NN=10000; % the test is achieved on 10000 OFDM symbols only. It is % possible to use all of the 100000 symbols, but it will % take more time. N=64; % number of subbands L=4; % oversampling factor C=16; % number of OFDM symbol candidates

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % phase factor matrix [B] generation p=[1 -1 j -j]; % phase factor possible values randn('state', 12345); B=randsrc(C,N,p); % generate N-point phase factors for each one of the % C OFDM candidates %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

for i=1:NN % calculate papr of original ofdm time_domain_signal=abs(ifft([ofdm_symbol(i,1:32) zeros(1,(L-1)*N) ofdm_symbol(i,33:64)])); meano=mean(abs(time_domain_signal).^2); peako=max(abs(time_domain_signal).^2); papro(i)=10*log10(peako/meano);

    % B*ofdm symbol
    p=[];
    for k=1:C
        p=[p; B(k,:).*ofdm_symbol(i,:)];
    end

     % Transform Pi to Time Domain and find paprs
     for k=1:C
         pt(k,:)=abs(ifft([p(k,1:32) zeros(1,(L-1)*N) p(k,33:64)]));
         papr(k)=10*log10(max(abs(pt(k,:)).^2)/mean(abs(pt(k,:)).^2));
     end

    % find papr_min
    papr_min(i)=min(papr);
    end

figure [cy,cx]=ccdf(papro,0.1); semilogy(cx,cy) % CCDF of PAPR of the original OFDM hold on [cy,cx]=ccdf(papr_min,0.1); semilogy(cx,cy,'r') % CCDF of the modified OFDM (by SLM)

Answer 3

[Sombaran Gupta]

Can anyone explain me from the third line till the end of this code..?? a=imread('kmeans.jpg'); imshow(255-rgb2gray(a)); b=reshape(double(255-rgb2gray(a)),[],1); [idx,c]=kmeans(b,2,'emptyaction','singleton'); Fin_a=reshape(idx,[size(rgb2gray(a))]); [ind1, ind2]=find(Fin_a(:,:)==2); [ind3, ind4]=find(Fin_a(:,:)==1); figure(); imshow(255-rgb2gray(a)); hold on if(length(ind1)<length(ind3)) scatter(ind2,ind1,'g'); else scatter(ind4,ind3,'g');a end

Answer 4

[Tanzila Minhaj]

Can anyone help me understanding the code?

function [sigma] = constitutive_Armstrong_Frederick (e, E, sigma_y, H_k, H_nl)
% The function reconstructs the trend of the stress-strain curve
% according to the constitutive link of Armstrong & Frederick
% e = deformation vector
% E = Young's modulus [Pa]
% sigma_y = yield stress [Pa]
% H_k = instant kinematic work hardening parameter [Pa]
% H_nl = hardening parameter by Armstrong & Frederick [Pa]

% Data initialization
sigma (1) = 0;
alpha (1) = 0;
state (1) = 0;
for n = 1: length (e) -1
    % Stress calculation by elastic prediction
    sigma_e = sigma (n) + E * (s (n + 1) -e (n));
    % Increased effort assessed in the elastic prediction phase
    dsigma_e = sigma_e-sigma (n);
    % Purified elastic prediction of back stress
    sigma_e_tilde = sigma_e-alpha (n);
    % Evaluation of the yield function to less than a tolerance
    fi (n) = abs (sigma (n) -alpha (n)) - sigma_y;
    tol = 0.001;
    if abs (fi (n)) <toll
        fi (n) = 0;
    end
    % Recognition of the elastic / plastic phase by the function of
    % yield and the direction of deformation
    if fi (n) == 0 && (e (n + 1) -e (n)) * (e (n) -e (n-1))> 0
        % Plastic phase
        state (n + 1) = 1;
        % Sign of tensions (+ = traction, - = compression)
        sign = sign (dsigma_e);
        % Increase in the plastic dlambda multiplier
        dlambda = (a * sigma_e_tilde-sigma_y) / (a * (E + H_k) -H_nl * alpha (n));
        % Increase in plastic deformation
        de_pl = a * dlambda;
        % Increase in back stress
        Dalfa = H_k * de_pl-H_nl * abs (de_pl) * alpha (n);
        % Back stress
        alpha (n + 1) = alpha (n) + Dalfa;
        % Current voltage
        sigma (n + 1) = sign * sigma_y alpha + (n + 1);
    else
        % Elastic phase
        state (n + 1) = 0;
        sigma (n + 1) = sigma_e;
        de_pl = 0;
        alpha (n + 1) = alpha (n);
    end
end
end

Answer 5

[Ayshath Afra]

Can anyone help me understanding the code?

How to obtain the read cbir function

Answer 6

[Arun das H A]

clear all;
clf('reset');
cam = webcam(); % Create a webcam object
right = imread('RIGHT1.jpg');
left = imread('LEFT1.jpg');
noface = imread('no_face1.jpg');
straight = imread('STRAIGHT1.jpg');

detector = vision.CascadeObjectDetector(); % Create a face detector using Viola-Jones
detector1 = vision.CascadeObjectDetector('EyePairSmall'); % Create an eyepair detector

while true
    vid = snapshot(cam); % Get a snapshot from the webcam
    vid = rgb2gray(vid); % Convert to grayscale
    img = flip(vid, 2); % Flip the image horizontally
    
    bbox = step(detector, img); % Creating bounding boxes using detector
    
    if ~isempty(bbox)%if face exists
        biggest_box = 1;
        for i = 1:size(bbox, 1)%find the biggest face
            if bbox(i, 3) > bbox(biggest_box, 3)
                biggest_box = i;
            end
        end
        
        faceImage = imcrop(img, bbox(biggest_box, :)); % Extract the face from the image
        
        bboxeyes = step(detector1, faceImage); % locations of he eye pair using detector
        
        subplot(2,2,1), subimage(img); hold on; % Display full image
        for i = 1:size(bbox,1)%draw all the regions that contain face
            rectangle('position', bbox(i, :), 'lineWidth', 2, 'edgeColor', 'y');
        end
        
        subplot(2,2,3), subimage(faceImage); % Display face image
        
        if ~isempty(bboxeyes)%check if eyepair is available
            biggest_box_eyes = 1;
            for i = 1:size(bboxeyes, 1)%find the biggest eye pair
                if bboxeyes(i, 3) > bboxeyes(biggest_box_eyes, 3)
                    biggest_box_eyes = i;
                end
            end
            
            bboxeyeshalf = [bboxeyes(biggest_box_eyes, 1), bboxeyes(biggest_box_eyes, 2), bboxeyes(biggest_box_eyes, 3)/3, bboxeyes(biggest_box_eyes, 4)];%resize the eyepair width in half
            eyesImage = imcrop(faceImage, bboxeyeshalf(1, :)); % Extract the half eyepair from the face image
            
            eyesImage = imadjust(eyesImage); % Adjust contrast
            
            r = bboxeyeshalf(1, 4) / 4;
            [centers, radii, metric] = imfindcircles(eyesImage, [floor(r - r/4) floor(r + r/2)], 'ObjectPolarity', 'dark', 'Sensitivity', 0.93); % Hough Transform
            
            [M, I] = sort(radii, 'descend');
            eyesPositions = centers;
            
            subplot(2,2,2), subimage(eyesImage); hold on;
            viscircles(centers, radii,'EdgeColor','b');
            
            if ~isempty(centers)
                pupil_x = centers(1);
                disL = abs(0 - pupil_x);%disatance from left edge to center point
                disR = abs(bboxeyes(1,3)/3 - pupil_x);%distance right edge to center point
                subplot(2,2,4);
                if disL > disR + 16
                    subimage(right);
                elseif disR > disL
                    subimage(left);
                else
                    subimage(straight);
                end
            end
        end
    else
        subplot(2,2,4);
        subimage(noface);
    end
end

[Walter Roberson]

@Arun das H A

I do not understand how this code will help anyone to understand how the posted gaussfilt function works??