[SciPy-User] weave.inline version of scipy.stats.rankdata

[Keith Goodman]

scipy.stats.rankdata looks like an easy function to speed up with
weave.inline. Well, easy if I knew what I was doing. My attempt below
sometimes segfaults and sometimes gives the wrong result on the first
run (but gives the right result on subsequent runs). Can anyone spot
what I am doing wrong? (A note about the code: I used "if...else if"
since ||, or, wouldn't compile.)

import numpy as np
from scipy import weave
from scipy.weave import converters

def rankdata(a):

a = np.asarray(a, dtype=float)
n = len(a)
svec, ivec = fastsort(a)
newarray = np.zeros(n, float)

code = """
int dupcount;
float sumranks, averank;
sumranks = 0.0;
dupcount = 0;
for (int i=0; i < n; ++i) {
sumranks = sumranks + i;
dupcount = dupcount + 1;
if (i==n-1) {
averank = sumranks / float(dupcount) + 1;
for (int j=i-dupcount+1; j < i+2; ++j) {
newarray[ivec[j]] = averank;
}
sumranks = 0;
dupcount = 0;
}
else if (svec[i] != svec[i+1]) {
averank = sumranks / float(dupcount) + 1;
for (int j=i-dupcount+1; j < i+2; ++j) {
newarray[ivec[j]] = averank;
}
sumranks = 0;
dupcount = 0;
}
}
"""
err = weave.inline(code,
['svec', 'ivec', 'n', 'newarray'],
type_converters=converters.blitz,
compiler='gcc')

return newarray

def fastsort(a):
# fixme: the wording in the docstring is nonsense.
"""Sort an array and provide the argsort.

Parameters
----------
a : array

Returns
-------
(sorted array,
indices into the original array,
)
"""
it = np.argsort(a)
as_ = a[it]
return as_, it
_______________________________________________
SciPy-User mailing list
SciPy...@scipy.org
http://mail.scipy.org/mailman/listinfo/scipy-user