New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

include_ordered for tckgen and tckedit #1238

Closed
wants to merge 68 commits into
base: dev
from

Conversation

Projects
None yet
4 participants
@LeeReidCSIRO
Contributor

LeeReidCSIRO commented Jan 22, 2018

Hi Guys,

Why

For years I've been repeatedly frustrated by stray tracks looping when performing tractography.
For example, tracking
S1M1 --> Posterior Limb of the Internal Capsule --> Brainstem
one commonly came across tracks that went
S1M1 --> thalmus or other-cortex --> Brainstem --> Posterior Limb IC

There were always work arounds - additional exclusion ROIs, making assumptions about brain size and setting the max-length etc, but these were all very hack-ish and tended to be a bit of trouble when dealing with non-standard brains.

What

This pull request adds -include_ordered <image_path> to tckedit and tckgen. This parameter enforces these include regions to be passed through in a particular order.

For example, to track from A, to B, to C:
tckgen fod.mif -seed_image a.nii.gz -seed_unidirectional -include_ordered a.nii.gz -include_ordered b.nii.gz -include_ordered c.nii.gz track.tck

This can be combined with standard -include rois as needed. For example, to track from A, to B, to C and at any point pass through X:
tckgen fod.mif -seed_image a.nii.gz -seed_unidirectional -include_ordered a.nii.gz -include_ordered b.nii.gz -include_ordered c.nii.gz -include X.nii.gz track.tck

This argument requires the -seed_unidirectional flag to be set, because bidirectional tracking does not have defined start and end points. If this flag is not set, an exception is thrown (just as with other checks).

Testing

To see this in action, you can find a sample script with accompanying fa/fod/ROIs here:
https://cloudstor.aarnet.edu.au/plus/s/NJeE6wV0Gw5kWIg
You will need to edit the script point it to your pull/build

Unit Testing

While building this, I couldn't find any unit testing framework in the code to check that what I had written actually works (?). I didn't want to go assuming which framework you wanted (or I missed) so I've written a very basic one from scratch which tests that the classes I've added. Please note that is purposefully does not use assert() because it's only set up to run tests in series: usually unit tests will continue to run even if one fails, but assert() prevents this. It's also nice to be able to run tests without re-configuring to DEBUG and rebuilding the entire project... I've tried to make this expandable for future additional tests, but ultimately it's probably much smarter to move to an standard established framework that integrates with builds etc.
There is an executable under testing called testing_unit_tests_tractography which runs these tests.

Please note that existing tests for mrconvert, population_template and tckconvert do not pass but I suspect that has nothing to do with these changes.

The Code

My background is predominantly in C#. I ported mrtrix 2.12 tckgen to C# years ago for my PhD work, in which I implemented something very similar. For me, rewriting this functionality in C++ was like hammering in a nail with my forehead. Beyond the syntactical/ differences, C# has standardised code conventions that do not match those of this C++ codebase. I say this because, due to force of habit, I've probably broken a few of your conventions here and there - sorry about this. I tried to adhere to what I saw, but I've probably missed a few things. If you have a style enforcer (e.g. StyleCop) with your rules put into it, please let me know so I don't do this again!

There are probably warnings related to the code I've created but the g++ compiler on mingw spits out so many warnings like this:

In file included from ../core/mrtrix.h:35:0,
                 from ../core/cmdline_option.h:26,
                 from ../core/app.h:26,
                 from ../core/command.h:21,
                 from cmd/testing_diff_fixel.cpp:15:
../core/types.h:283:40: warning: optimization attribute on 'uint64_t std::abs(uint64_t)' follows definition but the attribute doesn't match [-Wattributes]
   FORCE_INLINE uint64_t abs (uint64_t x) { return x; }
                                        ^
In file included from C:/msys64/mingw64/include/c++/7.2.0/cstdlib:77:0,
                 from C:/msys64/mingw64/include/c++/7.2.0/stdlib.h:36,
                 from C:/msys64/mingw64/lib/gcc/x86_64-w64-mingw32/7.2.0/include/mm_malloc.h:27,
                 from C:/msys64/mingw64/lib/gcc/x86_64-w64-mingw32/7.2.0/include/xmmintrin.h:34,
                 from ../core/command.h:19,
                 from cmd/testing_diff_fixel.cpp:15:
C:/msys64/mingw64/include/c++/7.2.0/bits/std_abs.h:61:3: note: previous definition of 'long long int std::abs(long long int)' was here
   abs(long long __x) { return __builtin_llabs (__x); }
   ^~~

...that I have no real hope of finding those which belong to me! The fact that the unit tests pass and behaviour of the executables matches what I'd expect is enough for me to feel confident that this works as advertised.

I have documented my code to a moderate degree. A pull request I will put in in some months time (different feature) will have comments for some of the existing tractography code as well. The lack of commenting in the existing code made penetrating it quite difficult for someone who didn't write it. I hope the comments I've added help others in the future.

Lee

@thijsdhollander thijsdhollander requested review from jdtournier, thijsdhollander and Lestropie and removed request for thijsdhollander Jan 22, 2018

@thijsdhollander

This comment has been minimized.

Member

thijsdhollander commented Jan 22, 2018

Hi @LeeReidCSIRO, just a heads-up; completely apart from the contents of the pull request: when introducing features of the scale and impact of something like this, we typically create a pull request with the dev branch as basis. Such a pull request will also run the latest tests with the latest test data.

Please note that existing tests for mrconvert, population_template and tckconvert do not pass but I suspect that has nothing to do with these changes.

That would be odd... in principle changes can't make it into either master or dev unless tests pass. But we'll see what the results are once the actual tests complete via TravisCI of course.

@LeeReidCSIRO

This comment has been minimized.

Contributor

LeeReidCSIRO commented Jan 22, 2018

@thijsdhollander Ok, would you like me to try to merge with dev branch and resubmit this?

My assumption about the tests may be proven wrong... I wasn't aware that they should always build in order to in the master branch. I just had no recollection of any of these having to rebuild when I changed the code. If it's my fault I'll try to sort.

@Lestropie

This comment has been minimized.

Member

Lestropie commented Jan 22, 2018

Hi Lee,

I might not be able to go through this in detail for a little while (fellowship application 👎), but I can offer a few comments up front:

  • Is there actually a need to be able to define both ordered and non-ordered include regions in a single command execution? Because if not, a command-line interface more consistent with MRtrix3 may be to use the -include option as normal, but then have an additional command-line flag, e.g. -ordered_include, which enforces ordered visitation of the include ROIs.
    Note this approach would also simplify resolving this contribution with other updates I have on the way, which separates -include into -include_wm and -include_gm.

  • Rather than merging to MRtrix3:master, it will be necessary to instead rebase the pull request to MRtrix3:dev. Currently we only push changes to the master branch that are either overt bug fixes, or tagged updates that enhance capabilities / change behaviour. Since this is the latter, it will need to be first merged to dev (which is where all development code converges), and then it will form part of the following tag release.
    Edit: Beaten to it. You should be able to edit the PR by changing the target branch.

  • Don't stress too much about the code convention; we'll go through and clean this up as part of the PR. We've tried style enforcers in the past but never been content with them, so our standards are slightly ill-defined. If there are going to be more external contributions however me way revisit this.

  • Regarding the "optimisation attribute" warning: I tried to shut this up with 2a50847, but then reverted in f4a96a8. I've only seen these on G++7 on Windows, and am not quite sure what to make of them.

@thijsdhollander

This comment has been minimized.

Member

thijsdhollander commented Jan 22, 2018

Small update here: the tests (currently still those from master) do pass successfully on TravisCI. No worries there. @LeeReidCSIRO : I'll rebase the pull request for you to dev. You'll still need to merge dev into your "feature branch" of course, but the pull request will make that clear.

@thijsdhollander thijsdhollander changed the base branch from master to dev Jan 22, 2018

@thijsdhollander

This comment has been minimized.

Member

thijsdhollander commented Jan 22, 2018

There you go. You'll notice the pull request now telling you (at the bottom) that "This branch is out-of-date with the base branch", since the base branch is now dev.

@LeeReidCSIRO

This comment has been minimized.

Contributor

LeeReidCSIRO commented Jan 22, 2018

Thanks for all the help guys. Working on rebasing now...

@Lestropie in answer to your query:

Is there actually a need to be able to define both ordered and non-ordered include regions in a single command execution? Because if not, a command-line interface more consistent with MRtrix3 may be to use the -include option as normal, but then have an additional command-line flag, e.g. -ordered_include, which enforces ordered visitation of the include ROIs.

I see a few options here:

  1. include_gm/wm are never ordered:
    -include A -include_ordered B -include_gm C -include_wm D
    (this is simplest)
  2. include_gm/wm ordered flag applies to all:
    -include A -include B -include_gm C -include_wm D -ordered_include
  3. additional flag or attribute for each:
    -include A -include B -ordered -include_gm C -include_wm D -ordered
    or
    -include A -include B -ordered -include_gm C -include_wm -ordered D
    or
    -include A -include (B -ordered) -include_gm C -include_wm (-ordered D)
    or
    -include A -include B,ordered -include_gm C -include_wm ordered,D (I don't like use of commas etc like this: it's perfectly pareseable with a little work, but hard to read and open for misinterpretation if your file is named unfortunately)

Something like option 3 is my preference, because it gives the flexibility that it has now, and means you can add additional attributes when you need new features later on. e.g.
-include A -include B -ordered -include (-gm C) -include (-wm -ordered D)
I don't know what's pareasable though.

Is there ever a need to be able to do both ordered and unordered? I don't currently have a need, personally. I put things the way they were because behaviour is undefined if two ordered rois overlap. It is undefined because re-entry is not allowed.
I.E. for -include_ordered a -include_ordered b -include_ordered c
a-->b-->a-->c is illegal. Track rejected
a --> b --> a/b --> c ?? It has re-entered A, but it hasn't left b... NB the track will be rejected in the current implementation.
I can imagine awkward situations cropping up where a user cannot guarantee that ROIs do not share voxels (e.g. being moved from two different atlases or where anatomy genuinely overlaps), or a case where you'd like tracks to pass through X regions in order, while at some point passing through (say) the left of some midline or through delineated pathological tissue. As it stands, it can handle that situation nicely by setting one ROI to -include. Yes, you could do this in two steps - tckgen then tckedit - but this is a really annoying thing to do when you'd like a guaranteed number of streamlines and don't know the proportions that will be rejected. Is this a niche case? Probably.

I'm not sure what the include_gm / include_wm are were intended for so I can't judge how that comes into play regarding the options I put above. Thoughts?

Regarding implementation, if you have a peek at ROISet when you have time, you'll see that integrating include_gm / include_wm should be a piece of cake here. Essentially, these will either be added to ROISet (as an ROISubSet), or within ROISubSet itself, depending on whether these are ordered or otherwise.

@thijsdhollander

This comment has been minimized.

Member

thijsdhollander commented Jan 22, 2018

Ok, as always, I haven't fully thought this through (yet), but just throwing it out there before I forget:

You could make the ordering "constraint" on inclusion regions also available via an option -order that has a parameter that is an index list of the inclusion regions you want to be included in the order constraint, and the order they sit in. This confines the whole mechanism to a single option (because only particular users will need it in particular scenarios, so it'd be nice to confine it to that), but makes it even more flexible than the above, and very efficient at it. It would also mean that it doesn't matter where in the list of options the -order argument itself sits. It only relies on the order of the inclusion regions in the list of options as specified by the user. The index list parameter works as e.g. the second parameter to -coord in mrconvert, etc... you can essentially specify any list of indices in many different flexible formats.
Examples:

-include A -include B -include C -order 0,1,2 : constrains to A-B-C

-include A -include B -include C -order 0:2 : same thing

-include A -include B -include C -order 0:end : same thing

-include A -include B -include C -order 1,2,0 : constrains to B-C-A

-include A -include B -include C -order 1,2 : constrains to ...-B-...-C-... , so only B and C in order, and A can sit wherever it likes (but is of course still required to be included)

-include A -include B -include C -order end:-1:0 : constrains to C-B-A ; note how this allows you to easily track through a long list of include regions in exact order, and only have to modify a tiny bit to the -order option to test tracking in reverse through them

-include A -include B -include C -order 0,1,2,0 : constrains to A-B-C-A; i.e. the streamline has to go through A-B-C in order and subsequently return to A

Using that latter flexibility you can come up with about any complex ordering, including repetition or omission (from the order constraint) of any inclusion region, that you can imagine; and specify it in a short and efficient manner (that could be scripted in lots of cool ways).

Also note that, in all of the above, the -order option itself can sit wherever it likes, like any other option. All of this could also have any other option sitting in between without issues. The only "order" of options that matters here is that of the -include options relative to each other, and that only matters when using the -order option of course, where it has to be clearly defined where those indices "point" to.

Any of the -include options above could just as well become an -include_wm, -include_gm or who knows what other "type" of include we may come up with. In this way, the implementation of new "types" of -include regions is entirely disconnected (independent) from the ordering constraint mechanism (from the point of view of the user interface; but almost certainly equally so from the point of view of the implementation, although that is in principle less important for the end-user).

@thijsdhollander

This comment has been minimized.

Member

thijsdhollander commented Jan 22, 2018

because behaviour is undefined if two ordered rois overlap. It is undefined because re-entry is not allowed.

Just to comment on that in particular: it doesn't have to be like that, in principle? If re-entry is not allowed, you could still enter regions in order, right? E.g. if A and B overlap partially, you could still first enter A, and then further "down the track" enter B, either while still in A, or after exiting A.

@LeeReidCSIRO

This comment has been minimized.

Contributor

LeeReidCSIRO commented Jan 22, 2018

I think the -order 0:2 etc looks nice. If lines are dropped within a command it's very readable:

#not readable:
tckgen fod.mif.gz -seed_image seed.nii.gz -include start_point.nii.gz -include mid_point_wm.nii.gz -include midpoint_gm -include endpoint_gm -order 0,1,3 -gm 2,3 -wm 0,1 -exclude exclude_mask.nii.gz tracks.tck

#pretty nice
tckgen fod.mif.gz -seed_image seed.nii.gz \
-include start_point.nii.gz \
-include mid_point_wm.nii.gz \
-include midpoint_gm \
-include endpoint_gm \
-order 0,1,3 \
-gm 2,3 \
-wm 0,1 \
-exclude exclude_mask.nii.gz tracks.tck

Two concerns though. Firstly, if this syntax is to be applied to other input images things get confusing. e.g. in the following example:

tckgen fod.mif.gz -seed_image seed.nii.gz -include mid_point_wm.nii.gz -exclude start_point.nii.gz -ignore_header 1 -ordered 0 exclude_mask.nii.gz tracks.tck

which image does -ignore_header refer to? -ordered 0 refers to the third image in the command.... It's surely confusing to the user for the second image to be 1, with the third being 0. If they are numbered counting from the same arg though, the programming becomes a headache: you'd have to break encapsulation and check every other module to see how many strings in the argument list refer to things that can have a property attached (images/meshes/spheres/etc). An initial stage of argument parsing may be able to be changed to handle this somehow, and pass the information into each class using it?

Finally, a very minor thought: you may (??) find yourself inundated with "it doesn't work" messages because people started counting at 1, not 0, or vice versa.

because behaviour is undefined if two ordered rois overlap. It is undefined because re-entry is not allowed.

Just to comment on that in particular: it doesn't have to be like that, in principle? If re-entry is not allowed, you could still enter regions in order, right? E.g. if A and B overlap partially, you could still first enter A, and then further "down the track" enter B, either while still in A, or after exiting A.

You're right, in principle it could be like you say - so long as A is entered first it should be OK to be in A & B. The code may become a tad inefficient though - at the moment it's very lightweight. I'm not sure if this matters: Rob just recently has highlighted to me that performance of tckgen is a sensitive issue. I'm happy to rewrite how this works if your team is of the consensus that A-->A&B something that should be handled.

@LeeReidCSIRO

This comment has been minimized.

Contributor

LeeReidCSIRO commented Jan 22, 2018

LeeReidCSIRO and others added some commits Jan 23, 2018

Prevent erroneous command-line usage warning
In #1282, a terminal warning was introduced in an attempt to catch cases where users failed to provide enough input arguments to command-line options.
This was refined in #1293 to prevent erroneous warnings from being displayed when having an input argument begin with a dash character is entirely reasonable.
This change further refines the logic dictating when a terminal warning should be issued, by identifying the use of piped images with command-line option inputs and not issuing the warning in that circumstance.
Merge pull request #1342 from MRtrix3/piped_option_arguments_warning_fix
Prevent erroneous command-line usage warning
Documentation fixes / updates after 3.0_RC3 merge
- Fix erroneous link syntax in fixel directory format description.
- Based on colleague feedback, clarify distinction between TmpFileDir and ScriptTmpDir config file entries.
- As suggested in #1336, Clarify locale requirements for MRtrix image format header key-value fields.
Closes #1336.
mrview Connectome tool: Fix streamtube initialisation
Failure to initialise member bool Connectome::have_streamtubes could lead to memory issues, as reported in #1330.
Merge pull request #1343 from MRtrix3/docs_hotfix_3.0_RC3
Documentation fixes / updates after 3.0_RC3 merge
dwi2response: Fix for absent -mask and no gradient table
If a mask is not provided to the script via the -mask option, the script will automatically generate a brain mask using the dwi2mask command. However, this command was being invoked upon the input image as provided by the user. If this image did not contain a gradient table in the header, and the user was relying on the -grad or -fslgrad option to provide this information, the dwi2mask call would fail. This change ensures that the dwi2mask command is executed on the DWI that has been imported into the temporary directory, and therefore is guaranteed to contain a diffusion gradient table.

jdtournier and others added some commits Jun 8, 2018

Merge pull request #1362 from MRtrix3/connectomestats_stdev_fix
connectomestats: Fix stdev output
Wait on stochastic threaded loop
Take chaFollowing 905317e, apply same change to stochastic threaded loop: call wait() function on output of Thread::run(), in order to properly catch exceptions.
Merge pull request #1350 from MRtrix3/mrview_connectome_streamtubes
mrview Connectome tool: Fix streamtube initialisation
Merge pull request #1365 from MRtrix3/mrview_tractography_geometry_mu…
…ltiple

mrview: Allow mulitple uses of -tractography.geometry
run.command(): Fix app cleanup with exitOnError=False
In rare instances, run.command() may be run with the argument exitOnError=False, allowing the underlying command to fail without aborting the script. When this occurs, only a warning-level message is issued. However, due to misplacement of code, such a circumstance would lead to the temporary directory not being erased on script completion, even if the script goes on to completion. This fix ensures that forcing the temporary directory to be retained only occurs if the underlying command both fails and is not permitted to fail.
checK_syntax: Look for %zu in print calls
As reported in #1367, these identifiers can lead to segmentation faults in C-style print functions when running on Windows. Therefore any new code should be checked to ensure that it is never again used.
Merge pull request #1384 from MRtrix3/fixelcfestats_remove_warning
fixelcfestats: Remove erroneous warning
Merge pull request #1367 from MRtrix3/fix_zu_printf_format
MSYS2: fix segfault due to use of %zu printf format
Merge pull request #1388 from MRtrix3/run_command_cleanup_fix
run.command(): Fix app cleanup with exitOnError=False

@Lestropie Lestropie referenced this pull request Jul 25, 2018

Open

Sync mrview #1403

@Lestropie

This comment has been minimized.

Member

Lestropie commented Aug 16, 2018

@LeeReidCSIRO: GitHub is currently showing this as coming from "unknown repository". Has your repo somehow moved / become private? It might be necessary to make a new PR if something has changed at your end.

On reflection I'm not as fussed on the -ordered_include flag idea; clashes with future changes shouldn't really be that bad. Interface-wise, -include_ordered is probably most consistent both with what we have and with how such data are utilised in the code. We've tried to avoid having complex command-line option groupings with brackets etc. thus far, and so would probably rather not use such unless there's no other way around it.

One thing I do want to have a look at (preferably in my code editor, hence needing to check the PR source) is the ROISet class. Here you've expanded this to contain both ordered and unordered sets. However the potential for ordering can only meaningfully be applied in the case of include regions, not exclusions or masks. So I'm wondering if the functionality would better live as a ROIOrderedSet, with ordered include regions being the only instance of that class; this would also be more consistent with the command-line interface. But I'd rather play with this myself and see whether or not I think it can work / makes more sense, rather than making a naive modification request.

@LeeReidCSIRO

This comment has been minimized.

Contributor

LeeReidCSIRO commented Aug 22, 2018

@Lestropie. Massive face-palm. I had issues with my github account a month or so back and forgot that this pull request was still active. It got wiped both on github and locally. Github has still kept all the diffs, though... If we want changes, and you don't have a local copy, one option is for me to try to get these diffs from github and start a new pull request referencing this one. Another option is to pull it into a different branch of mrtrix and make edits there (though I lack rights). Sorry about this.

Naming: I am happy so long as you are 👍

Exclude/Seed regions:
Personally, I try to keep things simple/clean usually, and try to reuse things even if it means a bit of redundancy. If there's a measurable performance drop from this new layout then I'm sure it would be relatively easy to split up or add a few functions to skip the odd boolean check.
Might be worth comparing performance to know whether this is worth our time? It's possible there's a slight performance impact but I suspect that it would be negligible. Side note: I have found previously (with the current public MRtrix release) that using an exclude region can be remarkably expensive, even if it's a tiny file and the ROI is never hit. I've never investigated why though.

@LeeReidCSIRO

This comment has been minimized.

Contributor

LeeReidCSIRO commented Aug 23, 2018

@Lestropie:
Update

I've recovered the code but I can't link it back to this thread. https://github.com/LeeReidCSIRO/mrtrix3/tree/ordered_include. Would you like me to make this into a new pull request and close this one?

The example script public link expired. It can also now be found at https://cloudstor.aarnet.edu.au/plus/s/oingnHX61hcvZWC

@Lestropie

This comment has been minimized.

Member

Lestropie commented Aug 24, 2018

Massive face-palm.

All good; sorry it's taken me this long to get to it. Been hella hectic trying to extinguish spot fires and put together a workshop having lost 2 core devs.

Would you like me to make this into a new pull request and close this one?

That's probably best I think: will mean that the merge will appropriately reference your current fork, we can request / propose changes properly, etc.. Can still give a link to this PR in the description in order to link to prior discussions.

If there's a measurable performance drop from this new layout then I'm sure it would be relatively easy to split up or add a few functions to skip the odd boolean check.

I'm not too concerned about the relative performance of the two approaches, it's more the logic of the code layout. At some point I will have to resolve this contribution against my suite of changes including differentiation between GM and WM ROIs, special handling of GM ROIs in the case of ACT, mesh ROIs, probably some other stuff I've done years ago and forgotten about. So I'm trying to keep each new capability as modular and independent as possible. But I'll need to spend some dedicated time looking at the code properly to decide whether or not there's a good justification one way or the other. I'd be curious to know @jdtournier's opinion on the possible code layout since it's his code.

Side note: I have found previously (with the current public MRtrix release) that using an exclude region can be remarkably expensive, even if it's a tiny file and the ROI is never hit.

Worth adding as its own issue; pretty sure we have open access for any GitHub account to create issues on our repo, other non-developers have done so in the past. Worst case scenario: someone tests it, they don't find anything, and we close the issue. Best case scenario: there's a problem, it gets fixed, and you get credit for having raised it. That way it'll be better documented, and there's a reminder there for yourself or anyone else to explicitly test it.

@LeeReidCSIRO

This comment has been minimized.

Contributor

LeeReidCSIRO commented Aug 27, 2018

Discussion moved to #1441

@Lestropie

This comment has been minimized.

Member

Lestropie commented Aug 28, 2018

Closed by move to #1441.

@Lestropie Lestropie closed this Aug 28, 2018

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment