A better way to create a response function for CSD

This implementation includes the recursive calibration of the response function for CSD, Using an FA threshold is not a very robust method.It is dependent on the dataset (non-informed used subjectivity), and still depends on the diffusion tensor (FA and first eigenvector), which has low accuracy at high b-value. This function recursively calibrates the response function, for more information see Tax et al. NeuroImage 2014.

This implementation is slow only on one line (788) which calls real_sph_harm. I would appreciate any help to optimize this function. Thx!
I believe, but I might be wrong, that here you could do something like `B_dwi = real_sph_harm(0, n[m == 0], theta[:, None], phi[:, None])`. The coefficients associated with the other rows of the matrix are going to be set to zero anyways and the functions should be orthogonal, so the result should not change.

This would require changes in a few other places in the code, for example the CSD model might need accept the response function as N / 2 + 1 coefficients (N being the SH order) of a radially symmetric function. But that might be the right thing to do anyways.

Please let me know if I'm just wrong about this or if I'm not explaining it well.