Merge pull request #160 from KrisThielemans/SinglesQuadruples

Add singles and quadruples and improve naming for other event types

- Adds `singleEvents` and `quadrupleEvents` to `EventTimeBlock`
- Add related policies and `*AreStored` variables to `ScannerInformation`
- Add ScannerInformation.delayedEventsPolicy
- Remove "optional" from `delayedEvents` and `tripleEvents`
  (***WARNING***: breaks backward compatibility on C++ level)
- More consistent naming: (***WARNING***: breaks backward compatibility)
  - Rename promptCoincidencesAreStored to promptEventsAreStored, and same for delayeds
  - Renamed ScannerInformation.coincidencePolicy to promptEventPolicy
  - Rename CoincidencePolicy.rejectMultiples to rejectHigherMultiples

Author: KrisThielemans

cpp/helpers/petsird_analysis.cpp

@@ -173,13 +173,15 @@ main(int argc, char const* argv[])
               for (unsigned mtype1 = 0; mtype1 < num_module_types; ++mtype1)
                 {
                   const petsird::TypeOfModulePair mtype_pair{ mtype0, mtype1 };
+
+                  // This code would need work to be able to handle a list-mode file without prompts
                   const auto& prompt_events = event_time_block.prompt_events[mtype0][mtype1];
 
                   // count events
                   num_prompts += prompt_events.size();
-                  if (event_time_block.delayed_events)
+                  if (scanner.delayed_events_are_stored)
                     {
-                      num_delayeds += (*event_time_block.delayed_events)[mtype0][mtype1].size();
+                      num_delayeds += event_time_block.delayed_events[mtype0][mtype1].size();
                     }
 
                   if (print_events)

cpp/helpers/petsird_generator.cpp

@@ -243,9 +243,12 @@ get_scanner_info()
 
   set_detection_efficiencies(scanner_info);
 
-  scanner_info.coincidence_policy = petsird::CoincidencePolicy::kRejectMultiples;
-  scanner_info.delayed_coincidences_are_stored = false;
+  scanner_info.prompt_event_policy = petsird::CoincidencePolicy::kRejectHigherMultiples;
+  scanner_info.single_events_are_stored = false;
+  scanner_info.prompt_events_are_stored = true;
+  scanner_info.delayed_events_are_stored = false;
   scanner_info.triple_events_are_stored = false;
+  scanner_info.quadruple_events_are_stored = false;
 
   return scanner_info;
 }
@@ -349,9 +352,10 @@ main(int argc, char* argv[])
       petsird::EventTimeBlock time_block;
       time_block.time_interval.start = t * EVENT_TIME_BLOCK_DURATION;
       time_block.time_interval.stop = (t + 1) * EVENT_TIME_BLOCK_DURATION;
-      time_block.prompt_events.resize(1);
-      time_block.prompt_events[type_of_module].resize(1);
-      time_block.prompt_events[type_of_module][type_of_module] = prompts_this_block;
+      auto& prompt_events = time_block.prompt_events;
+      prompt_events.resize(1);
+      prompt_events[type_of_module].resize(1);
+      prompt_events[type_of_module][type_of_module] = prompts_this_block;
       writer.WriteTimeBlocks(time_block);
     }
   writer.EndTimeBlocks();

model/CoincidenceEvent.yml

@@ -0,0 +1,38 @@
+# All information about a single event specified as indices (i.e. discretized).
+SingleEvent: !record
+  fields:
+    # identifiers of detecting bin (see doc for DetectionBin)
+    detectionBin: DetectionBin
+    # arrival time w.r.t. start of the time-block (in ps) (see EventTimeBlock).
+    timeOffsetInTimeBlock:  uint
+
+# All information about a coincidence event specified as indices (i.e. discretized).
+CoincidenceEvent: !record
+  fields:
+    # identifiers of the two detecting bins (see doc for DetectionBin)
+    detectionBins: DetectionBin*2
+    # an index into the tofBinEdges field in the ScannerInformation
+    tofIdx: uint
+
+# All information about a triple coincidence event specified as identifiers (i.e. discretized).
+TripleEvent: !record
+  fields:
+    # identifiers of the three detecting bins (see doc for DetectionBin)
+    # Note that it is possible that 2 "detecting elements" (i.e. module/element in the ExpandedDetectionBin)
+    # are equal (e.g. for inter-crystal Compton events)
+    detectionBins: DetectionBin*3
+    # timing differences w.r.t. first event, stored as
+    # indices into the tofBinEdges field in the ScannerInformation
+    # Note: only 2, corresponding to the arrival time differences of the second and third detectionBins
+    # listed w.r.t. the first detectionBin
+    tofIndices: uint*2
+
+# All information about a quadruple coincidence event specified as indices (i.e. discretized).
+QuadrupleEvent: !record
+  fields:
+    # identifiers of the four detecting bins (see doc for DetectionBin)
+    detectionBins: DetectionBin*4
+    # timing differences w.r.t. first event, stored as
+    # indices into the tofBinEdges field in the ScannerInformation.
+    # Note: only 3, see documentation for TripleEvent
+    tofIndices: uint*3

model/Events.yml

@@ -15,21 +15,57 @@ Header: !record
 # types of timeBlocks
 TimeBlock: [EventTimeBlock, ExternalSignalTimeBlock, BedMovementTimeBlock, GantryMovementTimeBlock, DeadTimeTimeBlock, SinglesHistogramTimeBlock]
 
+# Type definition for a list of single events (for a module of a specific type)
+ListOfSingleEvents: SingleEvent*
+
 # Type definition for a list of coincidences, where all coincidences are between 2 modules of specific types
 ListOfCoincidenceEvents: CoincidenceEvent*
 
 # Type definition for a list of triples, where all triples are between 3 modules of specific types
 ListOfTripleEvents: TripleEvent*
 
+# Type definition for a list of triples, where all triples are between 3 modules of specific types
+ListOfQuadrupleEvents: TripleEvent*
+
+# time-block that stores the detected (gamma photon) events as lists.
+# An event time block optionally stores
+# - singles
+# - "prompt" coincidences (i.e. 2 photons detected within the coincidence window),
+# - "delayed"coincidences (i.e. 2 photons with arrival time within the coincidence window after subtracting a
+# (scanner dependent) offset for the second event)
+# - triples
+# - quadruples
+# Note: PETSIRD currently does not allow storing higher order coincidences. These are very unlikely to occur in practice.
+# It is expected that most scanners will only store prompt and delayed coincidences. If singles are stored, it is expected
+# that coincidences are not stored, but this is not guaranteed.
+#
+# See also  the "policy" fields in ScannerInformation.
+# For instance, if coincidences, triples and quadruples are stored, it is likely that ScannerInformation.promptEventPolicy = rejectHigherMultiples,
+# ScannerInformation.tripleEventPolicy = rejectHigherMultiples, and ScannerInformation.quadrupleEventPolicy = rejectHigherMultiples as well,
+# but this is not guaranteed.
 EventTimeBlock: !record
   fields:
    # time interval for this time block
    # We suggest that for each coincidence/multiple event, the earliest time of arrival
    # of a photon is used for deciding which time block the event is entered in.
    # However, the impact of this choice is negligible and therefore can be vendor dependent.
+   # WARNING: if singleEvents are recorded, the time interval currently has to be smaller
+   # than 2^32 picoseconds due to how time is stored, see SingleEvent.
    timeInterval: TimeInterval
-   # Nested vector with lists of prompts in this time block
-   # There is one separate list for every pair of module types.
+
+   # nested vector with lists of singles
+   # If ScannerInformation.singleEventsAreStored is False, the vector should have
+   # size 0. Otherwise, there is one separate list for every module type.
+   # To get the list of singles for a type of module where the `detectionBin`
+   # is in a module of the first type, use
+   #  singleEvents[typeOfModule]
+   # Constraint: (size(singleEvents) > 0) == ScannerInformation.singleEventsAreStored
+   # Constraint: (size(singleEvents) == 0) || (size(singleEvents) == ScannerInformation.numberOfModuleTypes)
+   singleEvents: ListOfSingleEvents*
+
+   # nested vector with lists of prompts in this time block
+   # If ScannerInformation.promptEventsAreStored is False, the vector should have
+   # size 0. If not, there is one separate list for every pair of module types.
    # To get the list of prompts for 2 types of modules where the first `detectionBin`
    # given in the `CoincidenceEvent` is in a module of the first type, use
    #  promptEvents[typeOfModule1][typeOfModule2]
@@ -38,24 +74,40 @@ EventTimeBlock: !record
    # but this is currently not enforced. Therefore, to check all coincidences between two
    # different types of modules, a user will have to check two lists, i.e. also
    #  promptEvents[typeOfModule2][typeOfModule1]
-   # Constraint: size(promptsEvents) == ScannerInformation.numberOfModuleTypes
-   # Constraint: size(promptsEvents[*]) == ScannerInformation.numberOfModuleTypes
+   # Constraint: (size(promptEvents) > 0) == ScannerInformation.promptEventsAreStored
+   # Constraint: (size(promptEvents) == 0) || (size(promptsEvents) == ScannerInformation.numberOfModuleTypes)
+   # Constraint: (size(promptEvents) == 0) || (size(promptsEvents[*]) == ScannerInformation.numberOfModuleTypes)
    promptEvents: ListOfCoincidenceEvents**
-   # optional nested vector with lists of delayed coincidences in this time block
+
+   # nested vector with lists of delayed coincidences in this time block
    # To get the list of delayeds for 2 types of modules, use
    #  delayedEvents[typeOfModule1][typeOfModule2]
    # See promptEvents for more information.
-   # Constraint: size(delayedEvents) == ScannerInformation.numberOfModuleTypes
-   # Constraint: size(delayedEvents[*]) == ScannerInformation.numberOfModuleTypes
-   delayedEvents: ListOfCoincidenceEvents**?
-   # optional nested vector with lists of triple coincidences in this time block
+   # Constraint: (size(delayedEvents) > 0) == ScannerInformation.delayedCoincidencesAreStored
+   # Constraint: (size(delayedEvents) == 0) || (size(delayedEvents) == ScannerInformation.numberOfModuleTypes)
+   # Constraint: (size(delayedEvents) == 0) || (size(delayedEvents[*]) == ScannerInformation.numberOfModuleTypes)
+   delayedEvents: ListOfCoincidenceEvents**
+
+   # nested vector with lists of triple coincidences in this time block
    # To get the list of triples for 3 types of modules, use
    #  tripleEvents[typeOfModule1][typeOfModule2][typeOfModule3]
    # See promptEvents for more information.
-   # Constraint: size(tripleEvents) == ScannerInformation.numberOfModuleTypes
-   # Constraint: size(tripleEvents[*]) == ScannerInformation.numberOfModuleTypes
-   # Constraint: size(tripleEvents[*][*]) == ScannerInformation.numberOfModuleTypes
-   tripleEvents: ListOfTripleEvents***?
+   # Constraint: (size(tripleEvents) > 0) == ScannerInformation.tripleEventsAreStored
+   # Constraint: (size(tripleEvents) == 0) || (size(tripleEvents) == ScannerInformation.numberOfModuleTypes)
+   # Constraint: (size(tripleEvents) == 0) || (size(tripleEvents[*]) == ScannerInformation.numberOfModuleTypes)
+   # Constraint: (size(tripleEvents) == 0) || (size(tripleEvents[*][*]) == ScannerInformation.numberOfModuleTypes)
+   tripleEvents: ListOfTripleEvents***
+
+   # nested vector with lists of quadruple coincidences in this time block
+   # To get the list of quadruples for 4 types of modules, use
+   #  quadrupleEvents[typeOfModule1][typeOfModule2][typeOfModule3][typeOfModule4]
+   # See promptEvents for more information.
+   # Constraint: (size(quadrupleEvents) > 0) == ScannerInformation.quadrupleEventsAreStored
+   # Constraint: (size(quadrupleEvents) == 0) || (size(quadrupleEvents) == ScannerInformation.numberOfModuleTypes)
+   # Constraint: (size(quadrupleEvents) == 0) || (size(quadrupleEvents[*]) == ScannerInformation.numberOfModuleTypes)
+   # Constraint: (size(quadrupleEvents) == 0) || (size(quadrupleEvents[*][*]) == ScannerInformation.numberOfModuleTypes)
+   # Constraint: (size(quadrupleEvents) == 0) || (size(quadrupleEvents[*][*][*]) == ScannerInformation.numberOfModuleTypes)
+   quadrupleEvents: ListOfQuadrupleEvents****
 
 ExternalSignalTypeEnum: !enum
   values:
@@ -116,7 +168,12 @@ DeadTimeTimeBlock: !record
 SinglesHistogram: uint64[singlesDetectionBin]
 
 # A time block that stores a singles histogram.
-# See SinglesHistogramLevelType for information.
+# Note that all detected singles (within the ScannerInformation.singlesHistogramEnergyBinEdges)
+# should be stored. This is in contrast to the EventTimeBlock.singleEvents, where optionally
+# only singles that are in coincidence are recorded (see SingleEventPolicy).
+#
+# See also SinglesHistogramLevelType for information.
+#
 # Constraint: if ScannerInformation.singlesHistogramLevel == SinglesHistogramLevelType.none, there should
 #     be no time blocks of this type in the stream.
 #     Otherwise, the union of SinglesHistogramTimeBlock.timeInterval should be at least as large as the union of all EventTimeBlocks.timeInterval

model/Protocol.yml

@@ -1,12 +1,23 @@
+# Type definition for how to encode how the scanner records single events (if it does).
+SingleEventPolicy: !enum
+  values:
+    # store all single events
+    - all
+
+    # only store those that are in coincidence with another single
+    - rejectIfNotInCoincidence
+
+    # other options, to be listed in the future
+    - other
 
 # Type definition for how to encode how the scanner handles multiple coincidences when recording the prompts.
 # Due to various effects (such as high count rate, prompt gammas), it is possible that multiple single
 # events are detected within the coincidence window. This type encodes some different ways
-# that this multiple events are handled, and recorded in the coincidence stream.
+# that this multiple events are handled, and recorded in the EventTimeBlock.promptEvents field.
 CoincidencePolicy: !enum
   values:
-    # multiples will be rejected
-    - rejectMultiples
+    # multiples (more than 2 events in the coincidence window) will be rejected (i.e., not stored in the EventTimeBlock.promptEvents field)
+    - rejectHigherMultiples
 
     # multiples will be stored as a sequence of all pairs, e.g. a triple leads to 3 pairs
     - multiplesAsAllCoincidences
@@ -18,6 +29,25 @@ CoincidencePolicy: !enum
     # other options, to be listed in the future
     - other
 
+
+# Type definition for how to encode how the scanner handles higher order coincidences when recording the triples.
+TripleEventPolicy: !enum
+  values:
+    # quadruples etc are not stored in the EventTimeBlock.tripleEvents field
+    - rejectHigherMultiples
+
+    # other options, to be listed in the future
+    - other
+
+# Type definition for how to encode how the scanner handles higher order coincidences when recording the quadruples.
+QuadrupleEventPolicy: !enum
+  values:
+    # quadruples etc are not stored in the EventTimeBlock.quadrupleEvents field
+    - rejectHigherMultiples
+
+    # other options, to be listed in the future
+    - other
+
 # Type definition for how single histograms are stored.
 # Many scanners are able to store singles accumulated in time intervals, i.e.
 # histogrammed (as opposed to a list of singles with their corresponding DetectionBin).
@@ -115,14 +145,39 @@ ScannerInformation: !record
     #             (singlesHistogramLevel == SinglesHistogramLevelType.none && size(singlesHistogramEnergyBinEdges) == 0)
     singlesHistogramEnergyBinEdges: BinEdges*
 
-    # Encode how the scanner handles multiple coincidences
-    coincidencePolicy: CoincidencePolicy
+    # Encode how the scanner handles single events in the EventTimeBlock.singleEvents field
+    singleEventPolicy: SingleEventPolicy
+
+    # Encode how the scanner handles multiple coincidences in the EventTimeBlock.promptEvents field
+    promptEventPolicy: CoincidencePolicy
 
-    # a flag to indicate of delayed coincidences are recorded in the stream
-    delayedCoincidencesAreStored: bool
+    # Encode how the scanner handles multiple coincidences in the EventTimeBlock.delayedEvents field
+    # It is expected that this will have value rejectHigherMultiples, but it is not guaranteed.
+    delayedEventPolicy: CoincidencePolicy
 
-    # a flag to indicate of triple events are recorded in the stream
+    # Encode how the scanner handles higher order coincidences in the EventTimeBlock.tripleEvents field
+    # This field is ignored if tripleEventsAreStored == false.
+    tripleEventPolicy: TripleEventPolicy
+
+    # Encode how the scanner handles higher order coincidences in the EventTimeBlock.quadrupleEvents field
+    # This field is ignored if quadrupleEventsAreStored == false.
+    quadrupleEventPolicy: QuadrupleEventPolicy
+
+    # a flag to indicate of single events are recorded in the stream (see EventTimeBlock)
+    # Note that this is separate from the singlesHistogramLevel, where singles are histogrammed in a time-block.
+    singleEventsAreStored: bool
+
+    # a flag to indicate of prompt coincidences are recorded in the stream (see EventTimeBlock)
+    promptEventsAreStored: bool
+
+    # a flag to indicate of delayed coincidences are recorded in the stream (see EventTimeBlock)
+    delayedEventsAreStored: bool
+
+    # a flag to indicate of triple events are recorded in the stream (see EventTimeBlock)
     tripleEventsAreStored: bool
 
+    # a flag to indicate of quadruple coincidences are recorded in the stream (see EventTimeBlock)
+    quadrupleEventsAreStored: bool
+
     # coincidence detection efficiencies
     detectionEfficiencies: DetectionEfficiencies

model/ScannerInformation.yml

@@ -94,7 +94,7 @@ def parserCreator():
                         prompt_events = time_block.value.prompt_events[mtype0][
                             mtype1]
                         num_prompts += len(prompt_events)
-                        if time_block.value.delayed_events is not None:
+                        if scanner.delayed_events_are_stored:
                             num_delayeds += len(
                                 time_block.value.delayed_events[mtype0]
                                 [mtype1])

model/TripleEvent.yml

@@ -204,9 +204,12 @@ def get_scanner_info() -> petsird.ScannerInformation:
     # Now added the efficiencies
     scanner.detection_efficiencies = get_detection_efficiencies(scanner)
 
-    scanner.coincidence_policy = petsird.CoincidencePolicy.REJECT_MULTIPLES
-    scanner.delayed_coincidences_are_stored = False
+    scanner.coincidence_policy = petsird.CoincidencePolicy.REJECT_HIGHER_MULTIPLES
+    scanner.single_events_are_stored = False
+    scanner.prompt_events_are_stored = True
+    scanner.delayed_events_are_stored = False
     scanner.triple_events_are_stored = False
+    scanner.quadruple_events_are_stored = False
 
     return scanner