oxDNA Export with Helix Group Pitch / Yaw Bug

[ajvetturini]

Hello,

There seems to be a bug where when you apply pitch / yaw to a helix group and export to oxDNA format, the cross-section of the structure gets directionally skewed. This is best seen in the image below where on the left there is no pitch / yaw and the honeycomb lattice looks as expected. However on the right is when 60 degrees of pitch was assigned to the helix group that there are two issues: 1) the 6HB loses the cross-sectional shape and the helices overlap one another and 2) each strand within the bundle is now starting at a different position in space (depth)

If you want to see the issue / see a Python script that generates the malformed structure let me know and I can share the script or setup a Zoom call to screen share!

[dave-doty]

Here's a minimal reproducible example:

import scadnano as sc

def create_design() -> sc.Design:
    helices = [sc.Helix(max_offset=100) for _ in range(2)]
    design = sc.Design(helices=helices, grid=sc.square)

    design.groups[sc.default_group_name].pitch = 45

    design.draw_strand(0, 0).to(21)
    design.draw_strand(1, 0).to(21)

    return design

if __name__ == '__main__':
    d = create_design()
    d.write_scadnano_file(directory='output_designs')
    d.write_oxdna_files(directory='oxdna')

This makes this scadnano file: oxdna_export_with_pitch.zip

Here is how the exported oxDNA files appear:

Whereas with pitch set to 0 instead of 45, here is how it appears:

Here they are superimposed on top of each other:

Here's another example, this time with length-42 strands instead of length-21:

This gives a clue to the problem: it appears that the center of the strand, rather than the 5' end, is being placed at the origin. (Actually even that's not totally accurate: the origin of the whole system seems to be midway between the centers of the two strands.) Furthermore, each individual helix is being rotated around that point, rather than the shared origin of the whole HelixGroup.

[dave-doty]

I'm a bit confused actually. If you inspect the .dat file for pitch=0, the z-coordinates start at 0 and go up to 7.7952571026062465, as one would expect for 0.332 nm/base * 20 bases / 0.8518 nm/oxDNA length unit. But oxView doesn't draw the strands as having their 5' ends at z-coordinate 0, assuming the intersection of the three colored axes is actually the origin.

So this seems like a combination of expected behaviors from both oxView and scadnano:

1. In oxView, a strand its 5' end at the origin seems to be translated (or the axes are somehow misleading and not showing the oxView origin)
2. scadnano's export is not rotating all helices in a HelixGroup around the same axis, but is rotating each around a different axis.

[dave-doty]

I see one source of the problem. In the function _oxdna_get_helix_vectors, the variable position is set as follows (where position on the right side is from the field Helix.position):

position = position + group.position

In other words, it is accounting for the fact that the group's position needs to be added to the position of the individual helix. However, it does not account for the fact that rotating the whole HelixGroup around an axis would also shift the position of any helix whose position was not at the origin.

For example, if a helix had position (x=0, y=1, z=0), i.e., it lies 1 nm below the origin in the views above (where positive y goes down), then rotating the whole HelixGroup 90 degrees clockwise means its new position should be (x=0, y=0, z=-1). (i.e., the vector pointing from the origin to the helix goes from pointing straight down to pointing straight left.)