Python Scripting -- In Blender!
Doing my best to work smart, and not hard...
Those of you that know me IRL, are aware that I work as a civil engineer in Canada -- typically on municipal type projects for sewers, pump stations, water treatment, and roadway design. That sort of stuff. Either on the ground, or under it (generally).
I've found myself recently working very very iteraviely on some concept designs in AutoDesks Civil3D (very powerful, modern day engineering tool) -- but it felt tedious, and clunky, and involved a fair amount of going back and forth w/ spreadsheets to review / tweak / confirm designs.
Enter: Blender, and Python!
Currently, I'm very much in the beginning phases of understanding Blenders python API that they've put together to allow devs / designers to create / manipulate objects -- and I gotta say that I'm pretty impressed with the quality of documentation that they've put together. Something that I could learn from for my own projects and scripting adventures -- same with @exhaust.
Here's the gist of my simple script that I hacked together as a bit of a proof of concept:
from mathutils import * from math import * import bpy D = bpy.data # Find Site Surface meshes = D.meshes mesh = 0 for object in meshes: if object.name == "Surface": mesh = object # Find Lots lots =  ROWs =  for poly in mesh.polygons: if poly.material_index == 1: lots.append(poly) else: ROWs.append(poly) print("Found a bunch of lots!") # print(lots) # Find lot edges shared w/ ROWs row_verts =  for face in ROWs: for vertex in face.vertices: if vertex not in row_verts: row_verts.append(vertex) print(row_verts) for lot in lots: print("Checking Lot No. " + str(lot.index)) counter = 0 row_pts =  lot_verts = list(lot.vertices) set_lot = set(lot_verts) set_row = set(row_verts) lot_ROW_pts = list((set_row-(set_row-set_lot))) print(lot_ROW_pts) if len(lot_ROW_pts) == 2: frontage_vector = (mesh.vertices[lot_ROW_pts].co - mesh.vertices[lot_ROW_pts].co) print("Frontage Vector:", frontage_vector) frontage_vector = frontage_vector.normalized() print("IP#1:",mesh.vertices[lot_ROW_pts].co) print("IP#2:",mesh.vertices[lot_ROW_pts].co) print("Frontage Unit Vector:", frontage_vector) san_service_location = Vector( ( (mesh.vertices[lot_ROW_pts].co.x+mesh.vertices[lot_ROW_pts].co.x)/2, (mesh.vertices[lot_ROW_pts].co.y+mesh.vertices[lot_ROW_pts].co.y)/2, (mesh.vertices[lot_ROW_pts].co.z+mesh.vertices[lot_ROW_pts].co.z)/2 ) ) print("SAN:",san_service_location) wat_service_location = san_service_location+frontage_vector*2 print("WAT:",wat_service_location) bpy.ops.mesh.primitive_cube_add(location=san_service_location) bpy.context.active_object.name = "SAN Service" bpy.ops.mesh.primitive_ico_sphere_add(location=wat_service_location) bpy.context.active_object.name = "WAT Service"
My "proof of concept" script above, generally works as follows:
- I've created a random "subdivision layout" consisting of private properties (shaded blue w/ a simple material), and public road right-of-ways (ROWs) (shaded grey w/ a simple material) with some undulating topography;
- The script loads up the mesh ("surface" for you C3D types), and iterates through all the available faces / polygons;
- If the polygon has a blue material (
material_index == 1), it adds it to the list of private properties;
- If the polygon has a grey material (
material_index == 0), it adds it to the list of public ROWs;
- This is obviously a huge simplification -- but I'm starting simple.
- If the polygon has a blue material (
- After that, it detects with vertices of each lot are shared w/ the public ROWs, so we can determine which side of the private property is the frontage -- this is where our utility services will come in;
- Currently, I've only implemented support for sites that have TWO common vertices (this would be typical of a mid-block property), and have yet to support any lots w/ 3 (corner lot) or 4 (lots w/ a rear easement / alley way);
- After that, we take the two vertices that define the lot frontage (which would, in the real world, have iron pins set into the ground), and create a normalized (unit)
- For those of you that don't remember your geometry / linear algebra (or whatever branch of math vectors are most used), a unit vector is simple a vector where the LENGTH = 1. That makes it easy for use to follow along the frontage, and move along in 1m, 2m, 0.4m, or any arbitrary increments.
- For starters, we set the "Sanitary service location" as being at the direct midpoint of the frontage, and set the "Water service location" 2m offset from that.
|MMCD S6||MMCD W2b|
|Standard detail used in BC depicting how sanitary service connections should be constructed||Standard detail used in BC depicting how water service connections should be constructed|
- Then we create some placeholder objects so we can see that it's working.
I'm pretty fired up that this seems to work nicely, and quickly.
Next few steps should be:
- Create gravity service pipes at ~1.2m below surface, and connect to collection / distribution main in ROW;
- Gravity sewers would connect at a minimum -1% grade;
- Water services can connect anywhere, as they utilize flexible tubing and don't flow by gravity;
- Create a gravity sewer collection system that connects to the sanitary service connection at every private property;
- Sewers should maintain a minimum depth of cover (frost protection);
- Add manholes at sewer intersections and / or every 120m in linear sewer;
- Disregard practicality of construction for now. I anticipate that some sewers will be really really deep;
- Coming up with different alternatives will come later;
- Create a water distribution system that connects to the water service connection at every private property;
- This one is easier, as it doesn't need to flow by gravity;
- Maintain depth of cover;
- Look for conflicts between pipes;
Long term goals:
- Implement tools for hydraulic design, rather than just "connecting the dots";
- Apply approximate sewage generation rates for each lot (80~100% of average water consumption of 360L/person/day);
- Determine required pipe diameter and slope to convey flows safely;
- Implement unit rates for cost estimation;
- Come up w/ various options for servicing:
- 100% gravity collection to global low point;
- Check for opportunities for pump stations if required;
- Determine if certain properties should be serviced via individual, privately owned, pump stations rather than a gravity connection (ex -- ocean-front properties that are lower than the adjacent ROW);
- Compatibility w/ LandXML schema for import / export to AutoDesk Civil3D.
Any other civil engineers out here on the chain? I'd be keen to hear your thoughts on what you think may be an effective way to automate 'preliminary design' tasks like this. I think it'll be VERY complicated to script something for detailed design -- but for the early stage "options discovery", I think there's quite a bit of potential.
And yes -- I know that the python code above is like very inefficient. If you've got some cool suggestions for how to improve, I'd be keen to hear it!