               int offset_storage = old_global_count + lane;
               for (int idx = lane; idx < total_warp; idx += Warp::STRIDE, offset_storage += Warp::STRIDE)
               {
+                if (offset_storage >= output_xyz.size) break;
                 float x = storage_X[storage_index + idx];
                 float y = storage_Y[storage_index + idx];
                 float z = storage_Z[storage_index + idx];

                storage_Z[storage_index + offset + l] = points[l].z;
              }

              PointType *pos = output_xyz.data + old_global_count + lane;
              for (int idx = lane; idx < total_warp; idx += Warp::STRIDE, pos += Warp::STRIDE)
              {
                float x = storage_X[storage_index + idx];
                float y = storage_Y[storage_index + idx];
                float z = storage_Z[storage_index + idx];
                store_point_type (x, y, z, pos);
              }

              bool full = (old_global_count + total_warp) >= output_xyz.size;

        }

        __device__ __forceinline__ void
        store_point_type (float x, float y, float z, float3* ptr, int offset) const
        {
          *(ptr + offset) = make_float3 (x, y, z);
        }
      };

