vector index of consecutive gap (NaN) lengths?

Hi
For a vector A with random, sometime consecutive gaps of NaN, I want to develop a vector B of same length A that will indicate the length of local consecutive gaps for every value in A. B would have zeros for non-NaN locations in A.
so for
A = [2 4 NaN 7 9 NaN NaN NaN 32 NaN NaN 8];
I'd get
B = [0 0 1 0 0 3 3 3 0 2 2 0];
Ideas? Speed is always a virtue.
Thanks!
Tom

2 comentarios

Ryan
Ryan el 27 de Jun. de 2012
do you have the image processing toolbox?
tom
tom el 27 de Jun. de 2012
yes.

Iniciar sesión para comentar.

 Respuesta aceptada

Sean de Wolski
Sean de Wolski el 27 de Jun. de 2012
And if you like Ryan's idea but don't like bwlabel because it's evil:
A = [2 4 NaN 7 9 NaN NaN NaN 32 NaN NaN 8];
CC= bwconncomp(isnan(A));
n = cellfun('prodofsize',CC.PixelIdxList);
b = zeros(size(A));
for ii = 1:CC.NumObjects
b(CC.PixelIdxList{ii}) = n(ii);
end

3 comentarios

Ryan
Ryan el 27 de Jun. de 2012
And without the evil function it is faster! I have never used bwconncomp (mostly because I am used to matrices and vectors and not cells), but this may be the start of it.
Sean de Wolski
Sean de Wolski el 27 de Jun. de 2012
BWCONNCOMP makes BWLABEL irrelevant for everything except LABEL2RGB! For that you have LABELMATRIX to convert from the output of BWCONNCOMP.
Anyway, yes, CC.PixelIdxList contains the indices you need to do most matrix manipulations easily and can be fed directly to REGIONPROPS, all while being faster!
tom
tom el 28 de Jun. de 2012
Thanks Sean, this works great, even with NaNs at edges. -Tom

Iniciar sesión para comentar.

Más respuestas (3)

Ryan
Ryan el 27 de Jun. de 2012
Thomas' answer is faster, but here is my go:
A = [2 4 NaN 7 9 NaN NaN NaN 32 NaN NaN 8];
idx = isnan(A);
[B n]= bwlabel(idx);
C = B;
prop = regionprops(idx,'Area');
area = cat(1,prop.Area);
for ii = 1:n
B(C == ii) = area(ii);
end
B
Andrei Bobrov
Andrei Bobrov el 28 de Jun. de 2012
Editada: Andrei Bobrov el 30 de Jun. de 2012
A = [2 4 NaN 7 9 NaN NaN NaN 32 NaN NaN 8];
a = isnan(A);
t1 = find([true, diff(a)~=0]);
N = diff(t1);
out = zeros(size(A));
V = regionprops(a,'PixelIdxList');
out(cat(1,V.PixelIdxList)) = cell2mat(arrayfun(@(x)x*ones(x,1),N(a(t1))','un',0));
OR
A = [2 4 NaN 7 9 NaN NaN NaN 32 NaN NaN 8];
a = isnan(A);
n1 = regionprops(a,'Area');
out = a + 0;
out(a) = cell2mat(arrayfun(@(x)x*ones(1,x),[n1.Area],'un',0));
ADD variant
a = isnan(A);
t = [true,diff(a)~=0];
k = diff(find([t,true]));
k2 = k.*a(t);
out = k2(cumsum(t));
Thomas
Thomas el 27 de Jun. de 2012
A very crude way.. pretty sure can be done better...
A = [2 4 NaN 7 9 NaN NaN NaN 32 NaN NaN 8];
A(~isnan(A))=0;
A(isnan(A))=1;
c=diff(A);
start=find(c==1)+1;
stop=find(c==-1)+1;
out=stop-start;
for ii=1:length(out)
A(start(ii):(stop(ii)-1))=out(ii);
end
A

5 comentarios

tom
tom el 28 de Jun. de 2012
Thanks. This works nicely unless NaNs are on the ends, which is something it needs to be able to handle. Ideas?
Thomas
Thomas el 28 de Jun. de 2012
Editada: Thomas el 28 de Jun. de 2012
Based on your request
A = [2 4 NaN 7 9 NaN NaN NaN 32 NaN NaN 8 NaN];
A(~isnan(A))=0;
A(isnan(A))=1;
c=diff(A);
start=find(c==1)+1;
stop=find(c==-1)+1;
if length(stop)<length(start)
stop=[stop start(end)+1]
end
out=stop-start;
for ii=1:length(out)
A(start(ii):(stop(ii)-1))=out(ii);
end
A
If you're going for speed here is how they stack up:
Ryan's Answer: 0.014344 sec
Andrei's Answer: 0.05369 sec
Sean's Answer: 0.003353 sec
Thomas' Answer: 0.000045 sec
tom
tom el 29 de Jun. de 2012
This still doesn't work if A(1)=NaN; Thanks for continuing to work on it!
another iteration here NaN can be first,last or anywhere in the middle.. 'hopefully'
A = [NaN 4 NaN 7 9 NaN NaN NaN 32 NaN NaN 8];
A(~isnan(A))=0;
A(isnan(A))=1;
c=diff(A);
start=find(c==1)+1;
stop=find(c==-1)+1;
if length(stop)<length(start)
stop=[stop start(end)+1];
end
if length(start)<length(stop)
start=[start(1)-1 start];
end
out=stop-start;
for ii=1:length(out)
A(start(ii):(stop(ii)-1))=out(ii);
end
A

Iniciar sesión para comentar.

Categorías

Más información sobre MATLAB en Centro de ayuda y File Exchange.

Etiquetas

Preguntada:

tom
el 27 de Jun. de 2012

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by