add python implementation of Suheil Inati's iterative coil sensitivity estimation

[grlee77]

I added a pure python implementation of Gadgetron's coil_map_2d_Inati_Iter , coil_map_3d_Inati_Iter as mentioned in gadgetron/gadgetron#222. I assume channels is the first array dimension as is done in the other functions within ismrmrd-python-tools.

The csm example was updated to also compute the maps via this method.

Note:
I used the numpy docstring convention as it is the most common in scientific python packages and I find it easier to read. This also causes common python IDEs such as spyder to format the help nicely:

My preference would be to update the other functions in ismrmrd-python-tools to also use this convention. I am willing to do this if there is not an objection, but can instead change my function to match the current package convention if that is preferred.

[hansenms]

@inati, would you like to review?

[inati]

@grlee77 I have a python implementation in my branch. It's somewhat different from what's in the gadgetron and is the basis for a paper (in progress). Please take a look here:
https://github.com/inati/ismrmrd-python-tools/blob/coil_combine/ismrmrdtools/coils.py

The logic of the algorithm is not very well explained in the code, but there's a neat new twist that makes it invariant to unitary transformations (those that keep the noise white) and much more robust to coil compression. It also makes it very easy to parallelize. Once the paper is finished, I would like to change the gadgetron code to match this (both CPU and GPU versions). Not sure really how to proceed with this pull request right now. @hansenms and I need to discuss this a bit as well, but please take a look and let me know what you think. Perhaps we can also discuss on the call tomorrow.

[hansenms]

@inati, we can discuss this on the call tomorrow. However, if the code in the PR is correct, there is no need not to merge it, right? We can always replace it if and when a better implementation is available. There is no loss in merging code that does the right thing and looks OK?

[inati]

It looks to be a correct implementation of the algorithm from ISMRM2014 - a python port of a C++ port of the original matlab code - but I would need to check more carefully during daylight hours when my brain is working :-) Let's discuss a bit more before we merge and wait for @grlee77 to look at the other implementation.

[inati]

I'm going to merge this. Thank you @grlee77. I'm going to open an issue so that we can discuss how to improve our python docs/demos/webpage etc something like what the scikit folks are doing. I'm sure that @grlee77 will have some input there :-)

[grlee77]

Thanks @inati. I think the implementation in your branch is a bit easier to follow the intent with the various subfunctions separated out. The version in this PR has a couple nice features such as the threshold for stopping iterations and support for a non-isotropic smoothing kernel, that could also be easily incorporated there.

Does the version in your branch work for you on the example data that is in csm_estimation_demo.py? I called it via:

cs, rho, comim = calculate_csm_inati(coil_images)

but I got all NaN's in the outputs. Perhaps there is a divide by zero or something due to the noise-free data used in that demo. I will take a look.

[grlee77]

@inati:
as expected, my problem in the example was a divide by zero error due to the noise-free data. adding a small epsilon to the line computing csnorm fixes it:

i.e.

# normalize s*s*v, i.e. compute s*s
csnorm = np.sqrt(np.sum(np.abs(cs)**2,axis=0)) + eps
cs /= csnorm

[inati]

Correct! zero length vectors bad.