Diamond Ore Generation attempt

int spheres_made = 0;

bool func_207803_a(int config_size, float x_start, float x_end, float y_start, float y_end, float z_start, float z_end, int xstart, int ystart, int zstart, int width, int length) {

    uint32_t state = config_size + x_start + x_end + y_start + y_end + z_start + z_end + xstart+ystart+zstart+width+length;

    // BitSets are probably terrible for such small data
    std::vector<bool> bitset(width * width * length);
    std::vector<float> afloat(config_size * 4);

    for (int k = 0; k < config_size; ++k)
    {
        float t = (float) k / (float) config_size;

        float half_radius = config_size / 16.F;
        float random_times_half_radius = rnd(&state) * half_radius;

        // Point XYZ along vein axis start to end
        // Veins are only on the XZ plane so why we bother with Y calculations at all...
        // I have no idea, it's always going to result in 0 change.
        afloat[k * 4 + 0] = lerp(t, x_start, x_end);
        afloat[k * 4 + 1] = lerp(t, y_start, y_end);
        afloat[k * 4 + 2] = lerp(t, z_start, z_end);

        float taper = sinf(3.14159F * t); // Smaller at ends, larger in middle
        float taper_min = taper + 1.0F;  // Add one so we taper between 1x to 2x
        float taper_times_radius = taper_min * random_times_half_radius; // Combined taper and random effect
        float radius_min = taper_times_radius + 1.F; // Don't want a 0 radius, minimum radius = 1?
        float radius = radius_min / 2.F; // Not sure why the division again (comensate for +1's maybe) but it keeps the size in check

        // Diameter of vein is maximum of config.size / 4 at this point.

        afloat[k * 4 + 3] = radius;
    }
    // Some kind of smoothing? Can't visually tell the difference at all
    for (int k = 0; k < config_size - 1; ++k)
    { // ^v^v For each non-elimated ring
        if (afloat[k * 4 + 3] > 0)
        {
            for (int k3 = k + 1; k3 < config_size; ++k3)
            {  // ^v^v For each non-elimated ring after, moving towards end
                if (afloat[k3 * 4 + 3] > 0)
                {
                    float xdiff = afloat[k * 4 + 0] - afloat[k3 * 4 + 0];
                    float ydiff = afloat[k * 4 + 1] - afloat[k3 * 4 + 1];
                    float zdiff = afloat[k * 4 + 2] - afloat[k3 * 4 + 2];

                    float rdiff = afloat[k * 4 + 3] - afloat[k3 * 4 + 3];

                    /// If difference in ring size is greater than separation distance,
                    // No idea what the significance of this is... best guess is that
                    // it tries to make ends more flat but I'm not noticing that
                    // Or maybe a slight optimisation to remove small rings covered by
                    // previous rings?
                    if (rdiff * rdiff > xdiff * xdiff + ydiff * ydiff + zdiff * zdiff)
                    {
                        // Eliminate smaller ring
                        if (rdiff > 0)
                            afloat[k3 * 4 + 3] = -1.0;
                        else
                            afloat[k * 4 + 3] = -1.0;
                    }
                }
            }
        }
    }

    spheres_made = 0;

    for (int k = 0; k < config_size; ++k)
    {
        float ring_radius = afloat[k * 4 + 3];

        if (ring_radius >= 0)
        {
            // XYZ Position at point k
            float pos_x = afloat[k * 4 + 0];
            float pos_y = afloat[k * 4 + 1];
            float pos_z = afloat[k * 4 + 2];

            // Sub radius from current position as a starting point.
            // The max check is useless, ring_radius should have never been
            // set greater than the width/2 of our bounding volume
            int x_first = max(floor(pos_x - ring_radius), xstart);
            int y_first = max(floor(pos_y - ring_radius), ystart);
            int z_first = max(floor(pos_z - ring_radius), zstart);

            // Get final position to iterate to...

            // However, ring_radius can't be negative here, so this position is always greater than xyz_first
            // Therefore this max call is useless? even if the xyz_first was meant to be xyz_start + width it
            // would still be wrong... I don't get this...

            int x_last = max(floor(pos_x + ring_radius), x_first);
            int y_last = max(floor(pos_y + ring_radius), y_first);
            int z_last = max(floor(pos_z + ring_radius), z_first);
            for (int x = x_first; x <= x_last; ++x)
            {
                float x_perc_ring_radius = (x + 0.5f - pos_x) / ring_radius;

                if (x_perc_ring_radius < 1) // Don't do anything if outside of ring
                {
                    for (int y = y_first; y <= y_last; ++y)
                    {
                        float y_perc_ring_radius = (y + 0.5f - pos_y) / ring_radius;

                        if (mag(x_perc_ring_radius, y_perc_ring_radius) < 1)
                        {
                            for (int z = z_first; z <= z_last; ++z)
                            {
                                float z_perc_ring_radius = (z + 0.5f - pos_z) / ring_radius;

                                if (mag(x_perc_ring_radius, y_perc_ring_radius, z_perc_ring_radius) < 1)
                                {
                                    // 3D array indexing formula (bitset is a 3d array flattened)
                                    int l2 = x - (xstart) + (y - ystart) * width + (z - zstart) * width * length;

                                    if (!bitset[l2])
                                    {
                                        bitset[l2] = true;
                                        sphere s{ V3F(x, y, z), .5F, 2 };
                                        g_mc_world.spheres.push_back(s);
                                        ++spheres_made;
                                    }
                                }
                            }
                        }
                    }
                }
            }
        }
    }

    return spheres_made > 0;
}

bool generate(V3F pos, int config_size, int test=6)
{
    // How it works summary:
    // Vein made from a position P, size K, and angle A. The radius will be K/8.
    // The vein spans a line segment contain K points beginning at P-(angle * radius)
    // to P+(angle*radius). For each of these points randomly generate a 'ring radius'.
    // Do some filtering (or not) on it. Lastly for each point, iterate the blocks in the
    // bounding box around it (dimensions of ring radius), checking if the distance
    // to that block is within ring radius units. Basically you iterate over a sphere for
    // each of the K points.

    uint32_t state = config_size + pos.x + pos.y + pos.z + test*test;

    float angle = rnd(&state) * 3.14159F;
    //angle = 0;
    float radius = config_size / 8.f;

    // No need to caclulate these twice...
    float sinAngle = sin(angle) * radius;
    float cosAngle = cos(angle) * radius;

    float x_start = pos.x + sinAngle;
    float x_end   = pos.x - sinAngle;
    float z_start = pos.z + cosAngle;
    float z_end   = pos.z - cosAngle;

    float y_start = pos.y - (rand() % 3);
    float y_end   = y_start;

    int half_radius = ceil(radius / 2 + 0.5);

    int x_bound_start = pos.x - ceil(radius) - half_radius;
    int y_bound_start = pos.y - 2            - half_radius;
    int z_bound_start = pos.z - ceil(radius) - half_radius;

    int bound_width = 2 * (ceil(radius) + half_radius);
    int bound_length = 4 + 2 * half_radius;

    for (int x = x_bound_start; x <= x_bound_start + bound_width; ++x)
        for (int z = z_bound_start; z <= z_bound_start + bound_width; ++z)
            return func_207803_a(config_size, x_start, x_end, y_start, y_end, z_start, z_end, x_bound_start, y_bound_start, z_bound_start, bound_width, bound_length);

    return false;
}