Add dead-time time-block

cpp/petsird_analysis.cpp

@@ -122,7 +122,7 @@ main(int argc, char const* argv[])
       if (std::holds_alternative<petsird::EventTimeBlock>(time_block))
         {
           auto& event_time_block = std::get<petsird::EventTimeBlock>(time_block);
-          last_time = event_time_block.start;
+          last_time = event_time_block.time_interval.stop;
           num_prompts += event_time_block.prompt_events.size();
           if (event_time_block.delayed_events)
             num_delayeds += event_time_block.delayed_events->size();

cpp/petsird_generator.cpp

@@ -310,7 +310,8 @@ main(int argc, char* argv[])
       const auto num_prompts_this_block = poisson(gen);
       const auto prompts_this_block = get_events(header, num_prompts_this_block);
       petsird::EventTimeBlock time_block;
-      time_block.start = t * header.scanner.event_time_block_duration;
+      time_block.time_interval.start = t * header.scanner.event_time_block_duration;
+      time_block.time_interval.stop = (t + 1) * header.scanner.event_time_block_duration;
       time_block.prompt_events = prompts_this_block;
       writer.WriteTimeBlocks(time_block);
     }

model/Protocol.yml

@@ -14,14 +14,13 @@ Header: !record
 
 # types of timeBlocks
 # TODO more types could be needed
-TimeBlock: [EventTimeBlock, ExternalSignalTimeBlock, BedMovementTimeBlock, GantryMovementTimeBlock]
+TimeBlock: [EventTimeBlock, ExternalSignalTimeBlock, BedMovementTimeBlock, GantryMovementTimeBlock, DeadTimeTimeBlock]
 
 EventTimeBlock: !record
   fields:
-   # start time since ExamInformation.startOfAcquisition in ms
-   # Note: duration is given by ScannerInformation.eventTimeBlockDuration
-   start: uint
-   # TODO encode end time?
+   # time interval for this time block
+   # Constraint: all event-time-blocks need to have the same duration, see ScannerInformation.even_time_block_duration
+   timeInterval: TimeInterval
    # list of prompts in this time block
    # TODO might be better to use !array
    promptEvents: CoincidenceEvent*
@@ -50,21 +49,58 @@ ExternalSignalType: !record
 
 ExternalSignalTimeBlock: !record
   fields:
-   # start time since ExamInformation.startOfAcquisition in ms
-   start: uint
+   # time interval for this time block
+   timeInterval: TimeInterval
    # refer to ExternalSignalType.id
    signalID: uint
    # Note for triggers, this field is to be ignored
    signalValues: float*
 
 BedMovementTimeBlock: !record
   fields:
-    # start time since ExamInformation.startOfAcquisition in ms
-    start: uint
+    # time interval for this time block
+    timeInterval: TimeInterval
     transform: RigidTransformation
 
 GantryMovementTimeBlock: !record
   fields:
-    # start time since ExamInformation.startOfAcquisition in ms
-    start: uint
+    # time interval for this time block
+    timeInterval: TimeInterval
     transform: RigidTransformation
+
+# A time-block that stores the dead-time information for the given interval
+# TODO require that the complete scanning interval is "covered"?
+# Dead-time is encoded as "alive_time_fraction", i.e. 1 - dead_time_fraction
+# A component-based model is used, where the total alive-time for an pair of DetectorIds
+# is computed as the product of the (singles) alive-time of each detector times the
+# coincidence-alive-time of the corresponding module-pair:
+# aliveTimeFraction(detectorId0, energyIdx0, detectorId1, energyIdx1) =
+#   singlesAliveTimeFraction[detectorId0, energyIdx0] *
+#   singlesAliveTimeFraction[detectorId1, energyIdx1] *
+#   moduleCoincidenceAliveTimeFraction[module(detectorId0), module(energyIdx0)]
+DeadTimeTimeBlock: !record
+  fields:
+    # time interval for this time block
+    timeInterval: TimeInterval
+    # Singles dead-time array for the time interval, stored as fractions of "alive time"
+    # (1 means no dead-time, 0 means no detected counts)
+    # If sizeOfEnergyIdxDimension == 1, the fraction is assumed to be the same
+    # for all energy windows.
+    # Constraint: sizeOfDetectorIdDimension == number of all possible unique DetectorIds
+    # Constraint: (sizeOfEnergyIdxDimension == 1) or (sizeOfEnergyIdxDimension == ScannerInformation.numberOfEnergyBins)
+    singlesAliveTimeFraction: !array
+      items: float
+      dimensions: [detectorId, energyIdx]
+    # coincidence dead-time array for 2 modules in coincidence, stored as fractions of "alive time"
+    # (1 means no dead-time, 0 means no detected coincidences)
+    # If the size of this 2D array is (1,1), it is assumed that the corresponding alive-fraction is the same for all modules.
+    # Constraint: size(moduleCoincidenceAliveTimeFraction, 0) == 1 or total number of modules of all types
+    # Constraint: size(moduleCoincidenceAliveTimeFraction, 1) == 1 or total number of modules of all types
+    # Constraint: this matrix has to be symmetric
+    moduleCoincidenceAliveTimeFraction: !array
+      items: float
+      dimensions: [moduleId0, moduleId1]
+
+  computedFields:
+    sizeOfDetectorIdDimension: size(singlesAliveTimeFraction, 'detectorId')
+    sizeOfEnergyIdxDimension: size(singlesAliveTimeFraction, 'energyIdx')

python/petsird_analysis.py

@@ -72,11 +72,11 @@ def parserCreator():
         last_time = 0
         for time_block in reader.read_time_blocks():
             if isinstance(time_block, petsird.TimeBlock.EventTimeBlock):
-                last_time = time_block.value.start
+                last_time = time_block.value.time_interval.stop
                 num_prompts += len(time_block.value.prompt_events)
                 if time_block.value.delayed_events is not None:
                     num_delayeds += len(time_block.value.delayed_events)
-                print("=====================  Events in time block from ",
+                print("=====================  Events in time block until ",
                       last_time, " ==============")
                 for event in time_block.value.prompt_events:
                     energy_1 += energy_mid_points[event.energy_indices[0]]

python/petsird_generator.py

@@ -257,12 +257,14 @@ def get_events(header: petsird.Header,
         header = get_header()
         writer.write_header(header)
         for t in range(NUMBER_OF_TIME_BLOCKS):
-            start = t * header.scanner.event_time_block_duration
+            time_interval = petsird.TimeInterval(
+                start=t * header.scanner.event_time_block_duration,
+                stop=(t + 1) * header.scanner.event_time_block_duration)
             num_prompts_this_block = rng.poisson(COUNT_RATE)
             prompts_this_block = list(
                 get_events(header, num_prompts_this_block))
             # Normally we'd write multiple blocks, but here we have just one,
             # so let's write a tuple with just one element
             writer.write_time_blocks((petsird.TimeBlock.EventTimeBlock(
-                petsird.EventTimeBlock(start=start,
+                petsird.EventTimeBlock(time_interval=time_interval,
                                        prompt_events=prompts_this_block)), ))