And here is the OpenSCAD version, credit to user3717023 for the answer:

/* [Hidden] Constants */
point_color="lime";
golden_ratio = (sqrt(5)+1)/2; //1.618

/* [Dimensions] */
boundary_diameter = 120;
point_diameter = 4;

/* [Spiral] */
point_count = 450;
alpha=0;
golden_angle_coeff = 1.0;//0.99995;// e.g. 0.9995;

boundary_radius = boundary_diameter/2;
points_radius = point_diameter/2;
golden_angle = (360 - (360 / golden_ratio)) * golden_angle_coeff; //aka 137.5;

function radius(k,n,b) = (k > n-b) ?  1 : sqrt(k-1/2)/sqrt(n-(b+1)/2);
    
module sunflower( point_count, outer_radius, angle_stride, alpha) 
{    
    boundary_points = round(alpha * sqrt(point_count));  
    for (k = [1:point_count])    {
        r = radius(k, point_count, boundary_points) * outer_radius;
        theta = k * angle_stride;          
        rotate([0,0,theta]) translate([r,0,0]) children();
    }
}

union() {
    cylinder(d=boundary_radius*2, h=3, $fn=50);
    sunflower(point_count, boundary_radius-points_radius, golden_angle, alpha)
        color(point_color) cylinder(r=points_radius, h=3.5);
}