Fieldmaps acquired during each part of a session splits into two not differentiated

[neurorepro]

It is relatively common to have a single theoretical session splits practically into two sessions due to the length of the protocol, so e.g. anat, task 1, task 2, fmap in the morning and DWI, task 3, fmap in the afternoon.

heudiconv can deal with this with the option -g all to group all files no matter the series UID which changes after each practical session and would give an error otherwise.

The issue I am having is that when the grouping is applied to a dataset with a fieldmap having the same series id (e.g. 2_myfieldmap) during each practical session, the two fieldmaps are not differentiated as can be seen in dicominfo.tsv. The resulting detected fieldmap magnitude ('M' in protocol type) has then twice the number of files, and same for the phase, resulting in this kind of output for conversion:

│       │   ├── sub-s01_ses-01_magnitude1.json
│       │   ├── sub-s01_ses-01_magnitude1.nii.gz
│       │   ├── sub-s01_ses-01_magnitude2.json
│       │   ├── sub-s01_ses-01_magnitude2.nii.gz
│       │   ├── sub-s01_ses-01_magnitude3.json
│       │   ├── sub-s01_ses-01_magnitude3.nii.gz
│       │   ├── sub-s01_ses-01_magnitude4.json
│       │   ├── sub-s01_ses-01_magnitude4.nii.gz
│       │   ├── sub-s01_ses-01_phasediff1.json
│       │   ├── sub-s01_ses-01_phasediff1.nii.gz
│       │   ├── sub-s01_ses-01_phasediff2.json
│       │   └── sub-s01_ses-01_phasediff2.nii.gz

Would there be any suggestion on how to tackle this issue ?

Would the only solution be to convert each practical session (i.e. with a unique study UID) separately and group them later in a single theoretical session ? This would be quite cumbersome and a shame considering the great -g option normally made for this kind of case.

[yarikoptic]

The issue I am having is that when the grouping is applied to a dataset with a fieldmap having the same series id (e.g. 2_myfieldmap) during each practical session, the two fieldmaps are not differentiated as can be seen in dicominfo.tsv.

d'oh, not nice indeed. Moreover I am surprised that you get files from multiple studies get into a single .nii.gz! smells like something to clear up with dcm2niix guru @neurolabusc . So, what happens when you run dcm2niix directly on those files -- I would have assumed that dcm2niix would convert them into separate ones. What magnet is that btw? (note: filed #625).

sidenote:

heudiconv can deal with this with the option -g all to group all files no matter the series UID which changes after each practical session and would give an error otherwise.

yes, there is also -g accession_number (by accession ID) and even -g custom, see https://github.com/nipy/heudiconv/blob/master/docs/heuristics.rst#grouping-string-or-groupingfiles-dcmfilter-seqinfo

[neurolabusc]

dcm2niix is working as intended. It seems like heudiconv needs to resolve the correct IntendedFor when using the -g override.

[neurorepro]

Yes indeed dcm2niix is working as intended (because it looks at study uid ?)

Accession id is the same for all sequences across the two study UIDs in my case, but in our case we do want both field maps in the same session, just not considered as a single fieldmap.

There are different approaches to this, one would be to enrich the dicominfo file with fields allowing to differentiate the fieldmaps. Then heuristics.py could rely on that using appropriate info keys. A study UID field in dicominfo would actually be enough in this case. Unintuitive to add it since it is supposed to be constant by default, but it becomes quite useful when grouping accros study UIDs.

[yarikoptic]

sorry for the noise -- I have misread ". The resulting detected fieldmap magnitude ('M' in protocol type) has then twice the number of files" as ending with "number of volumes", my bad.

I am still not clear though on how magnitude is different/treated differently from phasediff which seems just to get two actual volumes for the place of 1 filename (e.g. sub-s01_ses-01_phasediff1.nii.gz). As for 'enrich' - do you think something like #581 would work for you? (then I should really get back to that draft)

[neurorepro]

@yarikoptic just to be clear, both magnitude and phasediff are combining files from two different sequences (happening in two different sessions study UIDs), so both of them show in dicominfo as sequences with twice the number of files. But at time of conversion you are right they are separate volumes although with 4 magnitudes and 2 phases instead of two separate sequences with 2 magnitudes and 1 phase.

So in dicominfo one would expect 4 rows, two with M, and two with P with number of files n_M and n_P respectively. But what happens is that they are only two rows (as if there was a single fieldmap) with number of files 2n_M and 2n_P.

If the study UID was taken into account to detect individual sequences in dicominfo and related, then I believe this would resolve the problem. And it would make sense as a single sequence is never spread over two sessions (hate to say never but that's maybe one rare case it could be said i believe).

[neurorepro]

@yarikoptic yes indeed that the PR i saw some time ago which would probably help with this kind of issue !

[neurorepro]

@yarikoptic Now that i had time to go back to that issue, I found out this is because series are identified only from their instance number and protocol name.

So if two different studies acquired in the morning and afternoon are merged (--grouping all), two series with the same instance number and protocol name will be fusionned as a single series !

More precisely, this is because the series identification is not done with the DICOM field SeriesInstanceUID but with... a custom SeriesID object (!) defined as:

heudiconv/heudiconv/dicoms.py

Lines 159 to 168 in f4e4896

	class SeriesID(NamedTuple):
	series_number: int
	protocol_name: str
	file_studyUID: Optional[str] = None
	def __str__(self) -> str:
	s = f"{self.series_number}-{self.protocol_name}"
	if self.file_studyUID is not None:
	s += f"-{self.file_studyUID}"
	return s

First question is: why the DICOM Series UID was not used for that ? Couldn't we just replace this custom SeriesID object with simply the value of SeriesInstanceUID ? There is the StudyUID attached to that object i know, but maybe there is an easy way to uncouple it from that object ?

In any case meanwhile I made PR #684 to minimally change the code to include the series UID in the grouping key. This way two series with same instance number and protocol name belonging to two different studies UID are correctly disentangled.

heudiconv/heudiconv/dicoms.py

Lines 160 to 171 in c6c7106

	class SeriesID(NamedTuple):
	series_number: int
	protocol_name: str
	series_uid: str
	file_studyUID: Optional[str] = None
	def __str__(self) -> str:
	suid_hex = hashlib.md5(self.series_uid.encode('utf-8')).hexdigest()[:16]
	s = f"{self.series_number}-{self.protocol_name}-{suid_hex}"
	if self.file_studyUID is not None:
	s += f"-{self.file_studyUID}"
	return s

Note: to still make the series id string representation readable, i computed a hash truncated to the first 16 characters, but the full series UID could be used instead (may be less readable though)

[neurorepro]

@yarikoptic i just saw you merged the other PR, great ! Did you also see #684 ? I am tagging you again as this issue is old so not sure you got the notification.

[yarikoptic]

I saw it but yet to provide feedback