| 1117 | } // anonymous namespace |
| 1118 | |
| 1119 | Expected<VoxelBitSet> segmentVolumeByGraphCut( const SimpleVolume & densityVolume, float k, const VoxelBitSet & sourceSeeds, const VoxelBitSet & sinkSeeds, ProgressCallback cb ) |
| 1120 | { |
| 1121 | MR_TIMER; |
| 1122 | |
| 1123 | if ( !reportProgress( cb, 0.0f ) ) |
| 1124 | return unexpectedOperationCanceled(); |
| 1125 | |
| 1126 | VoxelGraphCut vgc( densityVolume.dims ); |
| 1127 | vgc.resize( sourceSeeds, sinkSeeds ); |
| 1128 | if ( !reportProgress( cb, 2.0f / 16 ) ) |
| 1129 | return unexpectedOperationCanceled(); |
| 1130 | |
| 1131 | vgc.setupCapacities( densityVolume, k, sourceSeeds, sinkSeeds ); |
| 1132 | if ( !reportProgress( cb, 3.0f / 16 ) ) |
| 1133 | return unexpectedOperationCanceled(); |
| 1134 | |
| 1135 | vgc.setupNeighbors(); |
| 1136 | if ( !reportProgress( cb, 4.0f / 16 ) ) |
| 1137 | return unexpectedOperationCanceled(); |
| 1138 | |
| 1139 | auto parts = vgc.getSubtasks(); |
| 1140 | // parallel threads shall be able to safely modify elements in bit-sets |
| 1141 | const auto voxelsPerPart = ( int( vgc.getFullSpan().end ) / ( numSubtasks * BitSet::bits_per_block ) ) * BitSet::bits_per_block; |
| 1142 | |
| 1143 | auto sp = subprogress( cb, 4.0f / 16, 1.0f ); |
| 1144 | for ( int p = 0; p <= power; ++p ) |
| 1145 | { |
| 1146 | ParallelFor( parts, [&] ( size_t i ) |
| 1147 | { |
| 1148 | auto & part = parts[i]; |
| 1149 | if ( parts.size() == numSubtasks ) |
| 1150 | { |
| 1151 | part.span.begin = SeqVoxelId( i * voxelsPerPart ); |
| 1152 | part.span.end = ( i + 1 ) < numSubtasks ? SeqVoxelId( ( i + 1 ) * voxelsPerPart ) : vgc.getFullSpan().end; |
| 1153 | } |
| 1154 | |
| 1155 | VoxelGraphCut::Context context |
| 1156 | { |
| 1157 | .span = part.span, |
| 1158 | .stat = part.stat |
| 1159 | }; |
| 1160 | if ( parts.size() > 1 ) |
| 1161 | vgc.cutOutOfSpanNeiNeighbors( context ); |
| 1162 | vgc.buildForest( context, parts.size() == numSubtasks ); |
| 1163 | (void)vgc.segment( context ); // cannot fail without callback |
| 1164 | if ( parts.size() > 1 ) |
| 1165 | vgc.restoreCutNeighbor( context ); |
| 1166 | part.stat = context.stat; |
| 1167 | //part.stat.log( fmt::format( " after [{}, {})", part.span.begin, part.span.end ) ); |
| 1168 | } ); |
| 1169 | if ( parts.size() <= 1 ) |
| 1170 | { |
| 1171 | parts[0].stat.log( " final" ); |
| 1172 | break; |
| 1173 | } |
| 1174 | VoxelGraphCut::Statistics total; |
| 1175 | for ( size_t i = 0; 2 * i + 1 < parts.size(); ++i ) |
| 1176 | { |
no test coverage detected