Matlab - Entender funcion de Matlab

Hola! tengo que entender esta funcion para poder pasarla a C++, alguien me podría ayudar a comprenderla??
Gracias, Un Saludo.
Patricia


function [L,S,T,maxes] = find_landmarks(D,SR,N)
% [L,S,T,maxes] = find_landmarks(D,SR,N)
% Make a set of spectral feature pair landmarks from some audio data
% D is an audio waveform at sampling rate SR
% L returns as a set of landmarks, as rows of a 4-column matrix
% {start-time-col start-freq-row end-freq-row delta-time}
% N is the target hashes-per-sec (approximately; default 5)
% S returns the filtered log-magnitude surface
% T returns the decaying threshold surface
% maxes returns a list of the actual time-frequency peaks extracted.
%
% REVISED VERSION FINDS PEAKS INCREMENTALLY IN TIME WITH DECAYING THRESHOLD
%
% 2008-12-13 Dan Ellis [email protected]

if nargin < 3; N = 7; end % 7 to get a_dec = 0.998

% The scheme relies on just a few landmarks being common to both
% query and reference items. The greater the density of landmarks,
% the more like this is to occur (meaning that shorter and noisier
% queries can be tolerated), but the greater the load on the
% database holding the hashes.
%
% The factors influencing the number of landmarks returned are:
% A. The number of local maxima found, which in turn depends on
% A.1 The spreading width applied to the masking skirt from each
% found peak (gaussian half-width in frequency bins). A
% larger value means fewer peaks found.
f_sd = 30;

% A.2 The decay rate of the masking skirt behind each peak
% (proportion per frame). A value closer to one means fewer
% peaks found.
%a_dec = 0.998;
a_dec = 1-0.01*(N/35);
% 0.999 -> 2.5
% 0.998 -> 5 hash/sec
% 0.997 -> 10 hash/sec
% 0.996 -> 14 hash/sec
% 0.995 -> 18
% 0.994 -> 22
% 0.993 -> 27
% 0.992 -> 30
% 0.991 -> 33
% 0.990 -> 37
% 0.98 -> 67
% 0.97 -> 97



% A.3 The maximum number of peaks allowed for each frame. In
% practice, this is rarely reached, since most peaks fall
% below the masking skirt
maxpksperframe = 5;

% A.4 The high-pass filter applied to the log-magnitude
% envelope, which is parameterized by the position of the
% single real pole. A pole close to +1.0 results in a
% relatively flat high-pass filter that just removes very
% slowly varying parts; a pole closer to -1.0 introduces
% increasingly extreme emphasis of rapid variations, which
% leads to more peaks initially.
hpf_pole = 0.98;

% B. The number of pairs made with each peak. All maxes within a
% "target region" following the seed max are made into pairs,
% so the larger this region is (in time and frequency), the
% more maxes there will be. The target region is defined by a
% freqency half-width (in bins)
targetdf = 31; % +/- 50 bins in freq (LIMITED TO -32..31 IN LANDMARK2HASH)

% .. and a time duration (maximum look ahead)
targetdt = 63; % (LIMITED TO <64 IN LANDMARK2HASH)

% The actual frequency and time differences are quantized and
% packed into the final hash; if they exceed the limited size
% described above, the hashes become irreversible (aliased);
% however, in most cases they still work (since they are
% handled the same way for query and reference).


verbose = 0;

% Convert D to a mono row-vector
[nr,nc] = size(D);
if nr > nc
D = D';
[nr,nc] = size(D);
end
if nr > 1
D = mean(D);
nr = 1;
end

% Resample to target sampling rate
targetSR = 8000;
if (SR ~= targetSR)
D = resample(D,targetSR,SR);
end

% Take spectral features
% We use a 64 ms window (512 point FFT) for good spectral resolution
fft_ms = 64;
fft_hop = 32;
nfft = round(targetSR/1000*fft_ms);
S = abs(specgram(D,nfft,targetSR,nfft,nfft-round(targetSR/1000*fft_hop)));
% convert to log domain, and emphasize onsets
Smax = max(S(:));
% Work on the log-magnitude surface
S = log(max(Smax/1e6,S));
% Make it zero-mean, so the start-up transients for the filter are
% minimized
S = S - mean(S(:));
% This is just a high pass filter, applied in the log-magnitude
% domain. It blocks slowly-varying terms (like an AGC), but also
% emphasizes onsets. Placing the pole closer to the unit circle
% (i.e. making the -.8 closer to -1) reduces the onset emphasis.
S = (filter([1 -1],[1 -hpf_pole],S')');


% Estimate for how many maxes we keep - < 30/sec (to preallocate array)
maxespersec = 30;

ddur = length(D)/targetSR;
nmaxkeep = round(maxespersec * ddur);
maxes = zeros(3,nmaxkeep);
nmaxes = 0;
maxix = 0;

%%%%%
%% find all the local prominent peaks, store as maxes(i,:) = [t,f];

%% overmasking factor? Currently none.
s_sup = 1.0;

% initial threshold envelope based on peaks in first 10 frames
sthresh = s_sup*spread(max(S(:,1:min(10,size(S,2))),[],2),f_sd)';

% T stores the actual decaying threshold, for debugging
T = 0*S;

for i = 1:size(S,2)-1
s_this = S(:,i);
sdiff = max(0,(s_this - sthresh))';
% find local maxima
sdiff = locmax(sdiff);
% (make sure last bin is never a local max since its index
% doesn't fit in 8 bits)
sdiff(end) = 0; % i.e. bin 257 from the sgram
% take up to 5 largest
[vv,xx] = sort(sdiff, 'descend');
% (keep only nonzero)
xx = xx(vv>0);
% store those peaks and update the decay envelope
nmaxthistime = 0;
for j = 1:length(xx)
p = xx(j);
if nmaxthistime < maxpksperframe
% Check to see if this peak is under our updated threshold
if s_this(p) > sthresh(p)
nmaxthistime = nmaxthistime + 1;
nmaxes = nmaxes + 1;
maxes(2,nmaxes) = p;
maxes(1,nmaxes) = i;
maxes(3,nmaxes) = s_this(p);
eww = exp(-0.5*(([1:length(sthresh)]'- p)/f_sd).^2);
sthresh = max(sthresh, s_this(p)*s_sup*eww);
end
end
end
T(:,i) = sthresh;
sthresh = a_dec*sthresh;
end

% Backwards pruning of maxes
maxes2 = [];
nmaxes2 = 0;
whichmax = nmaxes;
sthresh = s_sup*spread(S(:,end),f_sd)';
for i = (size(S,2)-1):-1:1
while whichmax > 0 && maxes(1,whichmax) == i
p = maxes(2,whichmax);
v = maxes(3,whichmax);
if v >= sthresh(p)
% keep this one
nmaxes2 = nmaxes2 + 1;
maxes2(:,nmaxes2) = [i;p];
eww = exp(-0.5*(([1:length(sthresh)]'- p)/f_sd).^2);
sthresh = max(sthresh, v*s_sup*eww);
end
whichmax = whichmax - 1;
end
sthresh = a_dec*sthresh;
end

maxes2 = fliplr(maxes2);

%% Pack the maxes into nearby pairs = landmarks

% Limit the number of pairs that we'll accept from each peak
maxpairsperpeak=3;

% Landmark is <starttime F1 endtime F2>
L = zeros(nmaxes2*maxpairsperpeak,4);

nlmarks = 0;

for i =1:nmaxes2
startt = maxes2(1,i);
F1 = maxes2(2,i);
maxt = startt + targetdt;
minf = F1 - targetdf;
maxf = F1 + targetdf;
matchmaxs = find((maxes2(1,:)>startt)&(maxes2(1,:)<maxt)&(maxes2(2,:)>minf)&(maxes2(2,:)<maxf));
if length(matchmaxs) > maxpairsperpeak
% limit the number of pairs we make; take first ones, as they
% will be closest in time
matchmaxs = matchmaxs(1:maxpairsperpeak);
end
for match = matchmaxs
nlmarks = nlmarks+1;
L(nlmarks,1) = startt;
L(nlmarks,2) = F1;
L(nlmarks,3) = maxes2(2,match); % frequency row
L(nlmarks,4) = maxes2(1,match)-startt; % time column difference
end
end

L = L(1:nlmarks,:);

if verbose
disp(['find_landmarks: ',num2str(length(D)/targetSR),' secs, ',...
num2str(size(S,2)),' cols, ', ...
num2str(nmaxes),' maxes, ', ...
num2str(nmaxes2),' bwd-pruned maxes, ', ...
num2str(nlmarks),' lmarks']);
end

% for debug return, return the pruned set of maxes
maxes = maxes2;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function Y = locmax(X)
% Y contains only the points in (vector) X which are local maxima

% Make X a row
X = X(:)';
nbr = [X,X(end)] >= [X(1),X];
% >= makes sure final bin is always zero
Y = X .* nbr(1:end-1) .* (1-nbr(2:end));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function Y = spread(X,E)
% Each point (maxima) in X is "spread" (convolved) with the
% profile E; Y is the pointwise max of all of these.
% If E is a scalar, it's the SD of a gaussian used as the
% spreading function (default 4).
% 2009-03-15 Dan Ellis [email protected]

if nargin < 2; E = 4; end

if length(E) == 1
W = 4*E;
E = exp(-0.5*[(-W:W)/E].^2);
end

X = locmax(X);
Y = 0*X;
lenx = length(X);
maxi = length(X) + length(E);
spos = 1+round((length(E)-1)/2);
for i = find(X>0)
EE = [zeros(1,i),E];
EE(maxi) = 0;
EE = EE(spos+(1:lenx));
Y = max(Y,X(i)*EE);
end

Hola Patricia.
Primeramente cuales son los datos de entrada: D, SR, N; necesito su valor numerico de estos para poder ejecutarlo.
LA LINEA: [L,S,T,maxes] = find_landmarks(D,SR,N)
Te explico que signfica esto:
1) Este es un archivo de matlab tipo function, por eso al iniciar el programa la primera palabra es function, la cual es una palabra reservada de matlab, por lo tanto es color azul.
2). El nombre de la funcion es find_landmarks, con este nombre se debe guardar el archivo.
3). Tiene 3 datos de entradas(mas conocidos como argumentos de entrada) las cuales son D,RS,N
4). Tiene 4 datos salida(resultados) por eso tiene que ir entre corchetes [L, S, T, maxes]
es decir los datos de salida son: L, S, T y maxes que mas adelante enel programa debe crearse.

5) Todas las lineas que empiezan en porcentaje, son comentarios para que el usuario mejor que hace el programa.

LA LINEA: if nargin < 3; N = 7; end
1). la funcion nargin guarda el número de argumentos de entrada. En este programa el numero de argumentos de entrada debe ser 3 porque hay 3 datos de entradas.
Esta linea lo que hace, es que si el numero de argumentos es menor que 3 entonces a la variable N se asigna el valor de 7.
y continua el programa.
..
.
.
.
.
Explicarte linea por linea seria muy tedioso. Te sugiero que estudies matlab básico, los manuales lo puedes descargar de interent para puedas entender el programa con mayor precision.

Saludos.
JOSE JEREMIAS CABALLERO
ASESOR DE PROYECTOS CON MATLAB
PROFESOR DE METODOS NUMERICOS CON MATLAB
PROGRAMADOR EN MATLAB
[email protected]

Muchas gracias por tu ayuda de veras, ya lo compile en Matlab y funciona perfectamente. La señal D tengo que usar una cancion de tipo .wav pero en Matlab lo he estado probando con un coseno mismo para saber lo que hace.
N he usado el valor de 7 y SR 8000.

He leido manuales pero me gustaria sobre todo saber lo que hacen estas lineas, que no lo entiendo muy bien, lo demas ya he llegado ha entenderlo,aver si me podrias hechar una mano:




%%%%%
%% find all the local prominent peaks, store as maxes(i,:) = [t,f];

%% overmasking factor? Currently none.
s_sup = 1.0;

% initial threshold envelope based on peaks in first 10 frames
sthresh = s_sup*spread(max(S(:,1:min(10,size(S,2))),[],2),f_sd)';

% T stores the actual decaying threshold, for debugging
T = 0*S;

for i = 1:size(S,2)-1
s_this = S(:,i);
sdiff = max(0,(s_this - sthresh))';
% find local maxima
sdiff = locmax(sdiff);
% (make sure last bin is never a local max since its index
% doesn't fit in 8 bits)
sdiff(end) = 0; % i.e. bin 257 from the sgram
% take up to 5 largest
[vv,xx] = sort(sdiff, 'descend');
% (keep only nonzero)
xx = xx(vv>0);
% store those peaks and update the decay envelope
nmaxthistime = 0;
for j = 1:length(xx)
p = xx(j);
if nmaxthistime < maxpksperframe
% Check to see if this peak is under our updated threshold
if s_this(p) > sthresh(p)
nmaxthistime = nmaxthistime + 1;
nmaxes = nmaxes + 1;
maxes(2,nmaxes) = p;
maxes(1,nmaxes) = i;
maxes(3,nmaxes) = s_this(p);
eww = exp(-0.5*(([1:length(sthresh)]'- p)/f_sd).^2);
sthresh = max(sthresh, s_this(p)*s_sup*eww);
end
end
end
T(:,i) = sthresh;
sthresh = a_dec*sthresh;
end

% Backwards pruning of maxes
maxes2 = [];
nmaxes2 = 0;
whichmax = nmaxes;
sthresh = s_sup*spread(S(:,end),f_sd)';
for i = (size(S,2)-1):-1:1
while whichmax > 0 && maxes(1,whichmax) == i
p = maxes(2,whichmax);
v = maxes(3,whichmax);
if v >= sthresh(p)
% keep this one
nmaxes2 = nmaxes2 + 1;
maxes2(:,nmaxes2) = [i;p];
eww = exp(-0.5*(([1:length(sthresh)]'- p)/f_sd).^2);
sthresh = max(sthresh, v*s_sup*eww);
end
whichmax = whichmax - 1;
end
sthresh = a_dec*sthresh;
end

maxes2 = fliplr(maxes2);

%% Pack the maxes into nearby pairs = landmarks

% Limit the number of pairs that we'll accept from each peak
maxpairsperpeak=3;

% Landmark is <starttime F1 endtime F2>
L = zeros(nmaxes2*maxpairsperpeak,4);

nlmarks = 0;

for i =1:nmaxes2
startt = maxes2(1,i);
F1 = maxes2(2,i);
maxt = startt + targetdt;
minf = F1 - targetdf;
maxf = F1 + targetdf;
matchmaxs = find((maxes2(1,:)>startt)&(maxes2(1,:)<maxt)&(maxes2(2,:)>minf)&(maxes2(2,:)<maxf));
if length(matchmaxs) > maxpairsperpeak
% limit the number of pairs we make; take first ones, as they
% will be closest in time
matchmaxs = matchmaxs(1:maxpairsperpeak);
end
for match = matchmaxs
nlmarks = nlmarks+1;
L(nlmarks,1) = startt;
L(nlmarks,2) = F1;
L(nlmarks,3) = maxes2(2,match); % frequency row
L(nlmarks,4) = maxes2(1,match)-startt; % time column difference
end
end

L = L(1:nlmarks,:);

if verbose
disp(['find_landmarks: ',num2str(length(D)/targetSR),' secs, ',...
num2str(size(S,2)),' cols, ', ...
num2str(nmaxes),' maxes, ', ...
num2str(nmaxes2),' bwd-pruned maxes, ', ...
num2str(nlmarks),' lmarks']);
end

% for debug return, return the pruned set of maxes
maxes = maxes2;