@@ -263,8 +263,31 @@ def trace_cubes(self, beam: Union[dict[str, Any], Beam]) -> list[sitk.Image]:
263263 for i , cube in enumerate (rad_depth_cubes ):
264264 rel_depths = lengths * rho [i ]
265265 rel_depths = np .cumsum (rel_depths , axis = 1 ) - rel_depths / 2.0
266- ix_assign = np .unravel_index (ix [ix_remember_from_tracing ], cube .shape , order = "F" )
267- cube [ix_assign ] = rel_depths [ix_remember_from_tracing ]
266+
267+ try :
268+ ix_assign = np .unravel_index (ix [ix_remember_from_tracing ], cube .shape , order = "F" )
269+ except (ValueError , IndexError ):
270+ logger .error (
271+ "Error in unraveling indices from raytracing. Trying to recover..." ,
272+ exc_info = True ,
273+ )
274+ tmp_ix = ix [ix_remember_from_tracing ]
275+ rel_depths = rel_depths [ix_remember_from_tracing ]
276+
277+ wrong_values = np .logical_or (tmp_ix < 0 , tmp_ix >= cube .size )
278+ tmp_ix = tmp_ix [~ wrong_values ]
279+ rel_depths = rel_depths [~ wrong_values ]
280+
281+ # Remove the wrong values
282+ ix_assign = np .unravel_index (tmp_ix , cube .shape , order = "F" )
283+ logger .info (
284+ "Recovered %d indices for radiological depth cube" ,
285+ np .count_nonzero (wrong_values ),
286+ )
287+ cube [ix_assign ] = rel_depths
288+ else :
289+ cube [ix_assign ] = rel_depths [ix_remember_from_tracing ]
290+
268291 rad_depth_cubes [i ] = sitk .GetImageFromArray (cube )
269292 rad_depth_cubes [i ].CopyInformation (self .cubes [i ])
270293
0 commit comments