The code looks good to me. A few general points about this: 

1. Do you want to have this integrated into ConstrainedSphericalDeconvModel? That is, have `response` be an optional input to the class and if it isn't defined, run this function and derive it. The downside of that is that it would make initialization of the Model class instance rather long, because of the fitting of the TensorModel here.

2. It is important to note that mrtrix takes a different approach. They do not fit a tensor, but instead fit the response function as a very low order (4?) CSD. It might not make much of a difference, but it will not give identical answers. 

3. Something to think about: It might be advantageous to focus on a particular part of the brain, instead of blindly taking all the high FA voxels in the volume. For example, integrating the automatic CC-finding algorithm to focus specifically on that part of the brain, where we know that there is one major fiber bundle passing through. This might help with the problems that people have been pointing out about this procedure.