added a simple program to export files in .vdb format#148
added a simple program to export files in .vdb format#148spyke7 wants to merge 23 commits intoMDAnalysis:masterfrom
Conversation
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## master #148 +/- ##
==========================================
+ Coverage 88.20% 88.23% +0.02%
==========================================
Files 5 6 +1
Lines 814 884 +70
Branches 107 122 +15
==========================================
+ Hits 718 780 +62
- Misses 56 62 +6
- Partials 40 42 +2 ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
|
@orbeckst , please review the OpenVDB.py file. After that, I will add some more test covering all the missing parts |
|
Thank you for contribution. I’m currently on holidays and will come back to reviewing open source contributions in the new year. Am 12/27/25 um 03:16 schrieb Shreejan Dolai ***@***.***>:spyke7 left a comment (MDAnalysis/GridDataFormats#148)
@orbeckst , please review the OpenVDB.py file. After that, I will add some more test covering all the missing parts
—Reply to this email directly, view it on GitHub, or unsubscribe.You are receiving this because you were mentioned.Message ID: ***@***.***>
|
orbeckst
left a comment
orbeckst
left a comment
There was a problem hiding this comment.
Thank you for your contribution. Before going further, can you please try your own code and demonstrate that it works? For instance, take some of the bundled test files such as 1jzv.ccp4 or nAChR_M2_water.plt, write it to OpenVDB, load it in blender, and show an image of the rendered density?
Once we know that it's working in principle, we'll need proper tests (you can look at PR #147 for good example of minimal testing for writing functionality).
CHANGELOG
Outdated
| Fixes | ||
|
|
||
| * Adding openVDB formats (Issue #141) |
There was a problem hiding this comment.
not a fix but an Enhancement – put it into the existing 1.1.0 section and add you name there.
There was a problem hiding this comment.
In the CHANGELOG, this PR and issue are in the 1.1.0 release, so should I add my name in the 1.1.0 release or remove those lines and put them in the new section?
There was a problem hiding this comment.
Yes, now move it to the new section above since we released 1.1.0.
gridData/OpenVDB.py
Outdated
| for i in range(self.grid.shape[0]): | ||
| for j in range(self.grid.shape[1]): | ||
| for k in range(self.grid.shape[2]): | ||
| value = float(self.grid[i, j, k]) | ||
| if abs(value) > threshold: | ||
| accessor.setValueOn((i, j, k), value) |
There was a problem hiding this comment.
This looks really slow — iterating over a grid explicitly. For a start, you can find all cells above a threshold with numpy operations (np.abs(g) > threshold) and then ideally use it in a vectorized form to set the accessor.
|
fixed the CHANGELOG and OpenVDB.py. I didn't get the time to work on the blender part due to exams. I will surely try do it! |
The first two are for naChR_M2_water.vdb and the last two are for 1jzv.vdb. |
|
Good that you're able to load something into Blender. From a first glance I don;t recognize what I'd expect but this may be dependent on how you render in Blender. As I already said on Discord: Try to establish yourself what "correct" means. Load the original data in a program where you can reliably look at it. ChimeraX is probably the best for looking at densities; it can definitely read DX. Btw, the M2 density should look similar to the blue "blobs" on the cover of https://sbcb.bioch.ox.ac.uk/users/oliver/download/Thesis/OB_thesis_2sided.pdf |
The first one is for 1jzv.vdb and second for nAChR_M2_water.vdb (not done the shading/coloring) Can I proceed with the tests part? |
|
Mentioned in the Discord but also bringing up here: In your current examples (most obvious with the pore) is that the axis is flipped so that X is "up" compared to atomic coordinates which would have Z as up. |
Thank you for the update! will try to fix this |
I think this fixes the axis.. |
|
Ideally we would see this alongside the atoms or density from MN as well - to double check alignment because you might also need to flip one of the X or Y axes. |
|
The scales might be different (larger or smaller by factors of 10) but you can just scale inside of Blender by that amount to align the scales, but we want to be double checking alignemnt and axes. |
|
I have first of all added the MolecularNode add-on as given in the https://github.com/BradyAJohnston/MolecularNodes, and imported the 1jzv.pdb. After that import the .vdb file and there was difference in size of two. So I made the size the .pdb bigger. The centers of both of them are same and I didn't flipped any of the axes in the ss provided. I wrote a small blender py script to compare bounding boxes of the pdb and vdb objects to verify centroids, extents and axis alignment- output - The centroids are almost same I guess... |
|
@spyke7 It's still not 100% clear from your screenshots - can you import with the pore instead as that is more clear? And when you are taking a screenshot it would be more helpful to have the imported density in the centre of the screen rather than mostly empty space. |
|
Looks like you are attempting a standalone export to |
|
If this functionality can be added directly to GDF then we can also take advantage of that in MN going forwards. |
Agreed. In addition to exporting to |
This is a good point and something to consider as well. As far as I am aware Blender / MN (and other 3D animation packages) might be the only ones who use If there is anything out there that does take |
|
Yes, sort of: you need to add additional explicit keyword arguments to the top-level |
|
See #149 (comment) for a discussion for why we want to have explicit keywords. |
|
@PardhavMaradani can you please explain the center variable, more? As I cannot understand what to do for this.. |
The Here is an example of a density file (apbs.dx.gz) that is centered (left) and not (right): 
Here is a front view of the above: 
Here is a snippet from the code I pointed out in a previous comment: if center:
offset = -np.array(grid.shape) * 0.5 * gobj.delta
else:
offset = np.array(gobj.origin)
# apply transformations
vdb_grid.transform.preScale(np.array(gobj.delta) * world_scale)
vdb_grid.transform.postTranslate(offset * world_scale)If |
|
Thinking about this a bit more - @BradyAJohnston , given that the centering and scaling are just world transforms, do we really need to impose this upon |
orbeckst
left a comment
orbeckst
left a comment
There was a problem hiding this comment.
Minor changes, please run black on the files to get all formatting consistent.
Regarding the transformations it actually looks reasonable to me, but I want to hear more from @BradyAJohnston and @PardhavMaradani .
gridData/tests/test_vdb.py
Outdated
| assert tmpdir.join("auto.vdb").exists() | ||
|
|
||
| def test_write_vdb_with_metadata(self, tmpdir): | ||
| data = np.ones((3, 3, 3), dtype=np.float32) |
There was a problem hiding this comment.
Could use grid345 and then add metadata.
gridData/tests/test_vdb.py
Outdated
|
|
||
| class TestVDBWrite: | ||
| def test_write_vdb_from_grid(self, tmpdir, grid345): | ||
| data,g = grid345 |
gridData/tests/test_vdb.py
Outdated
| got = acc.getValue((i, j, k)) | ||
| assert got == pytest.approx(float(data[i, j, k])) | ||
|
|
||
| def test_write_vdb_default_grid_name(self, tmpdir): |
gridData/tests/test_vdb.py
Outdated
|
|
||
| voxel_size = grid_vdb.transform.voxelSize() | ||
|
|
||
| spacing=[voxel_size[0], voxel_size[1], voxel_size[2]] |
gridData/tests/test_vdb.py
Outdated
| vdb_field.write(outfile) | ||
|
|
||
| grids, metadata = vdb.readAll(outfile) | ||
| assert grids[0].name == 'direct_test' |
gridData/tests/test_vdb.py
Outdated
|
|
||
| spacing = [voxel_size[0], voxel_size[1], voxel_size[2]] | ||
|
|
||
| assert_allclose(spacing, [1.0, 2.0, 3.0], rtol=1e-5) |
There was a problem hiding this comment.
instead of [1.0, 2.0, 3.0] use the variable
| assert_allclose(spacing, [1.0, 2.0, 3.0], rtol=1e-5) | |
| assert_allclose(spacing, delta, rtol=1e-5) |
gridData/tests/test_vdb.py
Outdated
| assert acc.getValue((2, 3, 4)) == pytest.approx(5.0) | ||
| assert acc.getValue((7, 8, 9)) == pytest.approx(10.0) |
There was a problem hiding this comment.
instead of hard coding 5.0 and 10.0, access data
| assert acc.getValue((2, 3, 4)) == pytest.approx(5.0) | |
| assert acc.getValue((7, 8, 9)) == pytest.approx(10.0) | |
| assert acc.getValue((2, 3, 4)) == pytest.approx(data[2, 3, 4]) | |
| assert acc.getValue((7, 8, 9)) == pytest.approx(data[7, 8, 9]) | |
(and one could just make it a loop over index tuples if there were more than 2)
gridData/tests/test_vdb.py
Outdated
| grid_vdb = grids[0] | ||
| acc = grid_vdb.getAccessor() | ||
|
|
||
| assert acc.getValue((1, 1, 1)) == pytest.approx(1.0) |
gridData/OpenVDB.py
Outdated
| ] | ||
|
|
||
| vdb_grid.background = 0.0 | ||
| vdb_grid.transform = vdb.createLinearTransform(matrix) |
There was a problem hiding this comment.
I saw the tests (eg test_write_vdb_with_delta_matrix) that checked that reading the VDB file would reproduce the original delta and my understanding is that this works because of the transformations added here. I think it's quite important that we can roundtrip consistently so I would leave the transformations as they are as a default. (Correct me if I am wrong, please.)
If MN/Blender needs to scale/shift then we should make this possible on top of the default.
|
In this recent push, I have just applied the changes as asked in test_vdb.py. Will soon implement the scale and center in |
|
Yeah!, I added that 0.5 * delta offset because GDF uses a cell-centered convention. I will remove this, as this is just creating additional offset. |
There was a problem hiding this comment.
This looks all really good to me.
From my perspective, we only need to decide if the transformations should stay.
EDIT: Only just saw #148 (comment) — so we're keeping the transformation but remove the offset.
@BradyAJohnston @PardhavMaradani want some way to tweak the exports. Could you please leave a (blocking) review describing what you need to have added so that MN can make best use of the functionality?
|
Sorry should have time to look over this tomorrow. Adding the offset / centering on export is definitely something that could be handled by MN, but adding some transformation to the grid on export might still be useful more generally (or adding a transform as a Grid before export?). Will look over in more detail tomorrow. |
PardhavMaradani
left a comment
PardhavMaradani
left a comment
There was a problem hiding this comment.
Could you please leave a (blocking) review describing what you need to have added so that MN can make best use of the functionality?
Added what MN additionally needs (metadata support - apart from scale and center) and some general comments. Thanks
gridData/core.py
Outdated
| grid_name = self.metadata.get('name', 'density') | ||
|
|
||
| vdb_field = OpenVDB.OpenVDBField( |
There was a problem hiding this comment.
Do these have to be passed as params to __init__? (they are also currently marked as required params). This will have to be rewritten when import support is added as that will not have any of these values. I think it is better to keep the interfaces clean from the beginning. You can take a look at the mrc support on how this is handled for both cases. metadata will also need to be be available in the exported file.
There was a problem hiding this comment.
grid_name = self.metadata.get('name', 'density')
vdb_field = OpenVDB.OpenVDBField(
grid=self.grid,
origin=self.origin,
delta=self.delta,
name=grid_name
)The comment above is for the above lines...
gridData/OpenVDB.py
Outdated
|
|
||
| """ | ||
|
|
||
| def __init__(self, grid, origin, delta, name="density", tolerance=1e-10): |
There was a problem hiding this comment.
See previous comment about params here. About tolerance - shouldn't this be 0 by default so as to export the grid as is? How was this number determined and is this a generic value? My understanding is that this that OpenVDB sets any values around this value of the background to the background value. MN has a way to filter out the noise in a configurable way. I would try and avoid an arbitrary value as a default if possible. (If we know more why this was added - was this for reducing the file size, noise seen after import etc, we could see if there is a better solution)
There was a problem hiding this comment.
I have written for mainly reducing file size, or else, every value smaller than this will be there. though, it can be removed. It's definitely better to avoid an arbitrary value.
There was a problem hiding this comment.
Well. Should I remove the tolerance part completely from this?
as for BoolGrid, the function copyFromArray() does not accept tolerance, and thus inside prune(tolerance=False), should be present. So if I keep the tolerance part, then for float Grid, we need to mention it differently. Any thoughts?
There was a problem hiding this comment.
Set tolerance=None by default and treat None as the case where nothing is done to the data. If users want to change it then they can and if tolerance is not None or tolerance != 0 then run the pruning.
gridData/OpenVDB.py
Outdated
| if grid.ndim != 3: | ||
| raise ValueError(f"OpenVDB only supports 3D grids, got {grid.ndim}D") | ||
|
|
||
| self.grid = grid.astype(numpy.float32) |
There was a problem hiding this comment.
Is there a reason why everything is converted to float32's? Both GDF and OpenVDB support different grid data types, so we should make this generic? We will lose precision when grids have double values and be less memory efficient when we can use half grid etc.
gridData/OpenVDB.py
Outdated
|
|
||
| """ | ||
|
|
||
| vdb_grid = vdb.FloatGrid() |
There was a problem hiding this comment.
Similar to a comment above. OpenVDB supports different grid types and we should probably create one that corresponds to the grid data type?
There was a problem hiding this comment.
openvdb.GridTypes - gives this - [<class 'openvdb.FloatGrid'>, <class 'openvdb.BoolGrid'>, <class 'openvdb.Vec3SGrid'>],
So for this I guess - FloatGrid, BoolGrid, and Vec3SGrid these three are by deafult.
There was a problem hiding this comment.
From the same link above:
The Python module supports a fixed set of grid types. If the symbol PY_OPENVDB_WRAP_ALL_GRID_TYPES is defined at compile time, most of the grid types declared in openvdb.h are accessible in Python, otherwise only FloatGrid, BoolGrid and Vec3SGrid are accessible.
It looks like even the official module on conda-forge has only these fixed types:
>>> import openvdb
>>> openvdb.LIBRARY_VERSION
(13, 0, 0)
>>> openvdb.GridTypes
[<class 'openvdb.FloatGrid'>, <class 'openvdb.BoolGrid'>, <class 'openvdb.Vec3SGrid'>]
>>> hasattr(openvdb, "Int32Grid")
False
Blender packages its own version of openvdb and here is the output from Blender's Python Console:
>>> import openvdb
>>> openvdb.LIBRARY_VERSION
(12, 0, 0)
>>> openvdb.GridTypes
[<class 'openvdb.FloatGrid'>, <class 'openvdb.DoubleGrid'>, <class 'openvdb.BoolGrid'>, <class 'openvdb.Int32Grid'>, <class 'openvdb.Int64Grid'>, <class 'openvdb.Vec3SGrid'>, <class 'openvdb.Vec3IGrid'>, <class 'openvdb.Vec3DGrid'>, <class 'openvdb.PointDataGrid'>]
>>> hasattr(openvdb, "Int32Grid")
True
This is not ideal. MN has so far created the grids based on the corresponding data types (defaulting to float32 when there is no match) and this wasn't a problem because it runs within Blender. The hasattr checks are one way to check. There could be files with very different data types (the test nAChR_M2_water.plt file is float64 for example). I will defer to others on how best to deal with this. Thanks
There was a problem hiding this comment.
We'll have to work with what's available and this seems to be FloatGrid, BoolGrid (and Vec3SGrid, which we don't care about because all our densities are scalar).
Let's add a check that selects BoolGrid if the input array is a python (bool) or numpy bool numpy.bool and chooses FloatGrid for anything else. Add a note to the docs that limitations in OpenVDB can lead to loss of precision when input data is float64 (double) as the FloatGrid is float32 (single).
gridData/OpenVDB.py
Outdated
| # this is an explicit linear transform using per-axis voxel sizes | ||
| # world = diag(delta) * index + corner_origin | ||
| corner_origin = self.origin - 0.5 * self.delta |
There was a problem hiding this comment.
Since only the official openvdb module is being supported (I don't see pyopenvdb imports), it is probably simpler to use the preScale and postTranslate transforms directly.
There was a problem hiding this comment.
also, removed the corner_orgin, as the matrix is not needed anymore
There was a problem hiding this comment.
What's pyopenvdb – do you mean https://github.com/theNewFlesh/pyopenvdb ? It says it only supports Python 3.7 and 3.8 ... ???
There was a problem hiding this comment.
Inside of Blender openvdb was previously available as a module via import pyopenvdb but has been changed since 4.5 to import openvdb.
| vdb_grid.copyFromArray(self.grid, tolerance=self.tolerance) | ||
| vdb_grid.prune() | ||
|
|
||
| vdb.write(filename, grids=[vdb_grid]) |
There was a problem hiding this comment.
MN would need metadata support in the exported file - either all the current grid metadata or something explicitly passed during the export.
There was a problem hiding this comment.
please check the metadata part as implemented in the recent push
|
This might be a bit too late, but some thoughts on the design (maybe for the future) after reviewing the current implementation: The current |
|
Are you looking for a workflow such as the following @PardhavMaradani g = gdf.Grid("density.dx")
# make our VDB-like object that contains .vdb_grid as the VDB grid (eg FloatGrid)
gdf_vdb = gdf.OpenVDB.field(g)
# Then work with the VDB instance `vdb_grid`
gdf_vdb.vdb_grid.transform = createLinearTransform(matrix)
...If you provide code examples for how you would like to be able to use gdf then this would make things clear. |
Hi @orbeckst , yes, something along those lines. Here are some examples: Regular export from GDF: g = gdf.Grid("density.dx")
g.export("density.vdb", ...)The export options for above are any minimally required ones for basic functionality. Regular import into GDF: g = gdf.Grid("density.vdb")For someone like MN or others who want to add additional transforms and metadata to the vdb grid: g = gdf.Grid("density.dx")
vdb_grid = gdf.OpenVDB.grid_to_vdb(g)
vdb_grid.transform.preScale(...)
vdb_grid.transform.postTranslate(...)
vdb_grid["metadata_key_1"] = supported_type_value1
gdf.OpenVDB.write(vdb_grid, "/tmp/custom_grid.vdb", ...)
import openvdb
g1 = gdf.Grid("density1.dx")
g2 = gdf.Grid("density2.ccp4")
vdb_grid1 = gdf.OpenVDB.grid_to_vdb(g1)
vdb_grid2 = gdf.OpenVDB.grid_to_vdb(g2)
openvdb.write("/tmp/multiple_grids.vdb", grids=[vdb_grid1, vdb_grid2])The last two examples show how access to the openvdb grid can help with extensibility. Based on the above use cases,
The exporter could look something like: def _export_vdb(self, filename, ...):
...
gdf.OpenVDB.write(self, filename, ...)The importer could look like: def _load_vdb(self, filename, ...):
...
g = gdf.OpenVDB.read(filename, ...)
self._load(grid=g.grid, edges=g.edges, ...)Others who require additional |
gridData/OpenVDB.py
Outdated
| from gridData import OpenVDB | ||
| vdb_field = OpenVDB.field('density') | ||
| vdb_field.populate(grid, origin, delta) | ||
| vdb_field.write('output.vdb') |
gridData/OpenVDB.py
Outdated
| vdb_field = OpenVDB.field('density') | ||
| vdb_field.populate(grid, origin, delta) | ||
| vdb_field.write('output.vdb') | ||
|
|
|
Thanks for the use cases @PardhavMaradani, that's very helpful to see. We might be able to have gdf.OpenVDB contain simple functions and then introduce "convertors" for API interoperability (similar to what MDAnalysis is offering in the converter module. For instance, g = gdf.Grid("density.dx") # -> gdf.Grid
v = g.convert_to("vdb") # -> openvdb.GridBaseOnce we have this functionality, export is just doing this conversion before calling We can then also consider extending the converters to MRC objects. Eventually we could also add the functionality to drop OpenVDB or MRC objects into If we do a converter-style API then the gdf.OpenVDB module can be pretty light-weight because we don't really expect users to directly work with it. Does this sound like an interesting/clean way forward? |
|
@spyke7 I wanted to say that you're doing good work here! Don't be discouraged by the long discussions and the possibility that we want to change things again. You've demonstrated that the core of your code is working, now we can think about how this will best work long-term. Creating code that is actually used by people requires thought and discussion. The fact that we're having these discussions over your code means that this is something that we believe will have a long-term impact and is important enough to get right. |
yeah ofcourse, I will try the best updating the changes. And I can keep track of messages and reviews! Thanks. |
Using a generic |
Bringing this up from a review comment above to see if there is any possible way to address this as it seems a bit limiting.
|
4 participants
Hi @orbeckst
I have added
OpenVDB.pyinside gridData that simply export files in.vdbformat. Also I have addedtest_vdb.pyinside tests and it successfully passes.fix #141
Required Libraries -
openvdb
conda install -c conda-forge openvdbThere are many things that need to be updated like docs, etc, but I have just provided the file and test so that you can review it, and I can fix the problems. Please let me know if anything needs to be changed and updated.