1 #include "stepdata.h"
2 #include "atom_model.h"
3 #include "bond_model.h"
4 
5 using namespace Vipster;
6 
7 decltype (GUI::StepData::shader_map) GUI::StepData::shader_map;
8 
9 // TODO: update itself if applicable if scaling can be moved to shader
StepData(Step * step)10 GUI::StepData::StepData(Step* step)
11     : curStep{step}
12 {}
13 
StepData(StepData && dat)14 GUI::StepData::StepData(StepData&& dat)
15     : Data{std::move(dat)}
16 {
17     std::swap(atom_buffer, dat.atom_buffer);
18     std::swap(bond_buffer, dat.bond_buffer);
19     std::swap(cell_buffer, dat.cell_buffer);
20     std::swap(cell_mat, dat.cell_mat);
21     std::swap(curStep, dat.curStep);
22     std::swap(atRadFac, dat.atRadFac);
23     std::swap(object_map, dat.object_map);
24 }
25 
~StepData()26 GUI::StepData::~StepData()
27 {
28     for(auto &[context, objects]: object_map){
29         if(!objects.initialized) continue;
30         glDeleteVertexArrays(1, &objects.atom_vao);
31         glDeleteVertexArrays(1, &objects.bond_vao);
32         glDeleteVertexArrays(1, &objects.cell_vao);
33         glDeleteVertexArrays(1, &objects.sel_vao);
34         // create VBOs
35         glDeleteBuffers(1, &objects.atom_vbo);
36         glDeleteBuffers(1, &objects.bond_vbo);
37         glDeleteBuffers(1, &objects.cell_vbo);
38     }
39 }
40 
initGL(void * context)41 void GUI::StepData::initGL(void *context)
42 {
43     auto &shaders = shader_map[context];
44     if(!shaders.initialized){
45         auto &globals = global_map[context];
46         initAtomShader(globals, shaders);
47         initBondShader(globals, shaders);
48         initCellShader(globals, shaders);
49         initSelShader(globals, shaders);
50         shaders.initialized = true;
51     }
52     auto &objects = object_map[context];
53     if(!objects.initialized){
54         auto &globals = global_map[context];
55         // create VAOs
56         glGenVertexArrays(1, &objects.atom_vao);
57         glGenVertexArrays(1, &objects.bond_vao);
58         glGenVertexArrays(1, &objects.cell_vao);
59         glGenVertexArrays(1, &objects.sel_vao);
60         // create VBOs
61         glGenBuffers(1, &objects.atom_vbo);
62         glGenBuffers(1, &objects.bond_vbo);
63         glGenBuffers(1, &objects.cell_vbo);
64         // bind VBOs in VAO
65         initAtomVAO(globals, objects, shaders);
66         initBondVAO(globals, objects, shaders);
67         initCellVAO(globals, objects, shaders);
68         initSelVAO(globals, objects, shaders);
69         // mark as initialized
70         objects.initialized = true;
71     }
72     glBindVertexArray(0);
73 }
74 
initSelShader(GlobalContext & globals,ShaderContext & shaders)75 void GUI::StepData::initSelShader(GlobalContext &globals, ShaderContext &shaders)
76 {
77     auto &sel_shader = shaders.sel_shader;
78     sel_shader.program = loadShader(globals, "/select.vert", "/select.frag");
79     READATTRIB(sel_shader, vertex)
80     READATTRIB(sel_shader, position)
81     READATTRIB(sel_shader, vert_scale)
82     READATTRIB(sel_shader, hide)
83     READUNIFORM(sel_shader, pos_scale)
84     READUNIFORM(sel_shader, scale_fac)
85     READUNIFORM(sel_shader, offset)
86     READUNIFORM(sel_shader, pbc_instance)
87 }
88 
initSelVAO(GlobalContext & globals,ObjectContext & objects,ShaderContext & shaders)89 void GUI::StepData::initSelVAO(GlobalContext& globals, ObjectContext &objects, ShaderContext& shaders)
90 {
91     glBindVertexArray(objects.sel_vao);
92     auto &sel_shader = shaders.sel_shader;
93 
94     // sphere vertices
95     glBindBuffer(GL_ARRAY_BUFFER, globals.sphere_vbo);
96     glVertexAttribPointer(sel_shader.vertex, 3, GL_FLOAT, GL_FALSE, 0, nullptr);
97     glEnableVertexAttribArray(sel_shader.vertex);
98 
99     // atom positions
100     glBindBuffer(GL_ARRAY_BUFFER, objects.atom_vbo);
101     glVertexAttribPointer(sel_shader.position, 3,
102                           GL_FLOAT, GL_FALSE,
103                           sizeof(AtomProp), nullptr);
104     glVertexAttribDivisor(sel_shader.position, 1);
105     glEnableVertexAttribArray(sel_shader.position);
106 
107     // atom properties
108     glBindBuffer(GL_ARRAY_BUFFER, objects.atom_vbo);
109     glVertexAttribPointer(sel_shader.vert_scale, 1,
110                           GL_FLOAT, GL_FALSE,
111                           sizeof(AtomProp),
112                           reinterpret_cast<const GLvoid*>(offsetof(AtomProp, rad)));
113     glVertexAttribDivisor(sel_shader.vert_scale, 1);
114     glEnableVertexAttribArray(sel_shader.vert_scale);
115     glVertexAttribIPointer(sel_shader.hide, 1,
116                            GL_UNSIGNED_BYTE,
117                            sizeof(AtomProp),
118                            reinterpret_cast<const GLvoid*>(offsetof(AtomProp, hide)));
119     glVertexAttribDivisor(sel_shader.hide, 1);
120     glEnableVertexAttribArray(sel_shader.hide);
121 }
122 
initAtomShader(GlobalContext & globals,ShaderContext & shaders)123 void GUI::StepData::initAtomShader(GlobalContext& globals, ShaderContext &shaders)
124 {
125     auto &atom_shader = shaders.atom_shader;
126     atom_shader.program = loadShader(globals, "/atom.vert", "/atom.frag");
127     READATTRIB(atom_shader, vertex)
128     READATTRIB(atom_shader, position)
129     READATTRIB(atom_shader, vert_scale)
130     READATTRIB(atom_shader, color)
131     READATTRIB(atom_shader, hide)
132     READUNIFORM(atom_shader, offset)
133     READUNIFORM(atom_shader, pos_scale)
134     READUNIFORM(atom_shader, scale_fac)
135 }
136 
initAtomVAO(GlobalContext & globals,ObjectContext & objects,ShaderContext & shaders)137 void GUI::StepData::initAtomVAO(GlobalContext& globals, ObjectContext &objects, ShaderContext& shaders)
138 {
139     glBindVertexArray(objects.atom_vao);
140     auto &atom_shader = shaders.atom_shader;
141 
142     // sphere vertices
143     glBindBuffer(GL_ARRAY_BUFFER, globals.sphere_vbo);
144     glVertexAttribPointer(atom_shader.vertex, 3, GL_FLOAT, GL_FALSE, 0, nullptr);
145     glEnableVertexAttribArray(atom_shader.vertex);
146 
147     glBindBuffer(GL_ARRAY_BUFFER, objects.atom_vbo);
148     // atom positions
149     glVertexAttribPointer(atom_shader.position, 3,
150                           GL_FLOAT, GL_FALSE,
151                           sizeof(AtomProp),
152                           reinterpret_cast<const GLvoid*>(offsetof(AtomProp, pos)));
153     glVertexAttribDivisor(atom_shader.position, 1);
154     glEnableVertexAttribArray(atom_shader.position);
155 
156     // atom properties
157     glVertexAttribPointer(atom_shader.vert_scale, 1,
158                           GL_FLOAT, GL_FALSE,
159                           sizeof(AtomProp),
160                           reinterpret_cast<const GLvoid*>(offsetof(AtomProp, rad)));
161     glVertexAttribDivisor(atom_shader.vert_scale, 1);
162     glEnableVertexAttribArray(atom_shader.vert_scale);
163     glVertexAttribPointer(atom_shader.color, 4,
164                           GL_UNSIGNED_BYTE, GL_TRUE,
165                           sizeof(AtomProp),
166                           reinterpret_cast<const GLvoid*>(offsetof(AtomProp, col)));
167     glVertexAttribDivisor(atom_shader.color, 1);
168     glEnableVertexAttribArray(atom_shader.color);
169     glVertexAttribIPointer(atom_shader.hide, 1,
170                            GL_UNSIGNED_BYTE,
171                            sizeof(AtomProp),
172                            reinterpret_cast<const GLvoid*>(offsetof(AtomProp, hide)));
173     glVertexAttribDivisor(atom_shader.hide, 1);
174     glEnableVertexAttribArray(atom_shader.hide);
175 }
176 
177 
initBondShader(GlobalContext & globals,ShaderContext & shaders)178 void GUI::StepData::initBondShader(GlobalContext& globals, ShaderContext &shaders)
179 {
180     auto &bond_shader = shaders.bond_shader;
181     bond_shader.program = loadShader(globals, "/bond.vert", "/bond.frag");
182     READATTRIB(bond_shader, vertex)
183     READATTRIB(bond_shader, position)
184     READATTRIB(bond_shader, color1)
185     READATTRIB(bond_shader, color2)
186     READATTRIB(bond_shader, mMatrix)
187     READATTRIB(bond_shader, pbc_crit)
188     READUNIFORM(bond_shader, offset)
189     READUNIFORM(bond_shader, pos_scale)
190     READUNIFORM(bond_shader, pbc_cell)
191     READUNIFORM(bond_shader, mult)
192 }
193 
initBondVAO(GlobalContext & globals,ObjectContext & objects,ShaderContext & shaders)194 void GUI::StepData::initBondVAO(GlobalContext& globals, ObjectContext &objects, ShaderContext& shaders)
195 {
196     glBindVertexArray(objects.bond_vao);
197     auto &bond_shader = shaders.bond_shader;
198 
199     // cylinder vertices
200     glBindBuffer(GL_ARRAY_BUFFER, globals.cylinder_vbo);
201     glVertexAttribPointer(bond_shader.vertex,3,GL_FLOAT,GL_FALSE,0,nullptr);
202     glEnableVertexAttribArray(bond_shader.vertex);
203 
204     // model matrix (rotation)
205     glBindBuffer(GL_ARRAY_BUFFER, objects.bond_vbo);
206     glVertexAttribPointer(bond_shader.mMatrix, 3,
207                           GL_FLOAT,GL_FALSE,
208                           sizeof(BondProp), nullptr);
209     glVertexAttribDivisor(bond_shader.mMatrix, 1);
210     glVertexAttribPointer(bond_shader.mMatrix+1, 3,
211                           GL_FLOAT,GL_FALSE,
212                           sizeof(BondProp),
213                           reinterpret_cast<const GLvoid*>(offsetof(BondProp, mat[3])));
214     glVertexAttribDivisor(bond_shader.mMatrix+1, 1);
215     glVertexAttribPointer(bond_shader.mMatrix+2, 3,
216                           GL_FLOAT,GL_FALSE,
217                           sizeof(BondProp),
218                           reinterpret_cast<const GLvoid*>(offsetof(BondProp, mat[6])));
219     glVertexAttribDivisor(bond_shader.mMatrix+2, 1);
220     glEnableVertexAttribArray(bond_shader.mMatrix);
221     glEnableVertexAttribArray(bond_shader.mMatrix+1);
222     glEnableVertexAttribArray(bond_shader.mMatrix+2);
223 
224     // position
225     glVertexAttribPointer(bond_shader.position, 3,
226                           GL_FLOAT,GL_FALSE,
227                           sizeof(BondProp),
228                           reinterpret_cast<const GLvoid*>(offsetof(BondProp, pos)));
229     glVertexAttribDivisor(bond_shader.position, 1);
230     glEnableVertexAttribArray(bond_shader.position);
231 
232     // pbc conditions
233     glVertexAttribIPointer(bond_shader.pbc_crit, 4,
234                            GL_SHORT,
235                            sizeof(BondProp),
236                            reinterpret_cast<const GLvoid*>(offsetof(BondProp, mult)));
237     glVertexAttribDivisor(bond_shader.pbc_crit, 1);
238     glEnableVertexAttribArray(bond_shader.pbc_crit);
239 
240     // Color 1
241     glVertexAttribPointer(bond_shader.color1, 4,
242                           GL_UNSIGNED_BYTE,GL_TRUE,
243                           sizeof(BondProp),
244                           reinterpret_cast<const GLvoid*>(offsetof(BondProp, col_a)));
245     glVertexAttribDivisor(bond_shader.color1, 1);
246     glEnableVertexAttribArray(bond_shader.color1);
247 
248     // Color 2
249     glVertexAttribPointer(bond_shader.color2, 4,
250                           GL_UNSIGNED_BYTE,GL_TRUE,
251                           sizeof(BondProp),
252                           reinterpret_cast<const GLvoid*>(offsetof(BondProp, col_b)));
253     glVertexAttribDivisor(bond_shader.color2, 1);
254     glEnableVertexAttribArray(bond_shader.color2);
255 }
256 
initCellShader(GlobalContext & globals,ShaderContext & shaders)257 void GUI::StepData::initCellShader(GlobalContext& globals, ShaderContext &shaders)
258 {
259     auto &cell_shader = shaders.cell_shader;
260     cell_shader.program = loadShader(globals, "/cell.vert", "/cell.frag");
261     READATTRIB(cell_shader, vertex)
262     READUNIFORM(cell_shader, offset)
263 }
264 
initCellVAO(GlobalContext & globals,ObjectContext & objects,ShaderContext & shaders)265 void GUI::StepData::initCellVAO(GlobalContext& globals, ObjectContext &objects, ShaderContext& shaders)
266 {
267     glBindVertexArray(objects.cell_vao);
268 
269     auto &cell_shader = shaders.cell_shader;
270 
271     glBindBuffer(GL_ELEMENT_ARRAY_BUFFER, globals.cell_ibo);
272     glBindBuffer(GL_ARRAY_BUFFER, objects.cell_vbo);
273     glVertexAttribPointer(cell_shader.vertex, 3, GL_FLOAT, GL_FALSE, 0, nullptr);
274     glEnableVertexAttribArray(cell_shader.vertex);
275 }
276 
draw(const Vec & off,const PBCVec & mult,const Mat & cv,bool drawCell,void * context)277 void GUI::StepData::draw(const Vec &off, const PBCVec &mult,
278                          const Mat &cv, bool drawCell, void *context)
279 {
280     Vec tmp;
281     const auto &objects = object_map[context];
282     const auto &shaders = shader_map[context];
283     // atoms
284     glBindVertexArray(objects.atom_vao);
285     glUseProgram(shaders.atom_shader.program);
286     glUniform1f(shaders.atom_shader.scale_fac, atRadFac);
287     glUniformMatrix3fv(shaders.atom_shader.pos_scale, 1, 0, cell_mat.data());
288     for(int x=0;x<mult[0];++x){
289         for(int y=0;y<mult[1];++y){
290             for(int z=0;z<mult[2];++z){
291                 tmp = (off + x*cv[0] + y*cv[1] + z*cv[2]);
292                 glUniform3f(shaders.atom_shader.offset,
293                             static_cast<float>(tmp[0]),
294                             static_cast<float>(tmp[1]),
295                             static_cast<float>(tmp[2]));
296                 glDrawArraysInstanced(GL_TRIANGLES, 0,
297                                       atom_model_npoly,
298                                       static_cast<GLsizei>(atom_buffer.size()));
299             }
300         }
301     }
302     // bonds
303     glBindVertexArray(objects.bond_vao);
304     glUseProgram(shaders.bond_shader.program);
305     glUniform3i(shaders.bond_shader.mult, mult[0], mult[1], mult[2]);
306     glUniformMatrix3fv(shaders.bond_shader.pos_scale, 1, 0, cell_mat.data());
307     for(GLint x=0;x<mult[0];++x){
308         for(GLint y=0;y<mult[1];++y){
309             for(GLint z=0;z<mult[2];++z){
310                 tmp = (off + x*cv[0] + y*cv[1] + z*cv[2]);
311                 glUniform3f(shaders.bond_shader.offset,
312                             static_cast<float>(tmp[0]),
313                             static_cast<float>(tmp[1]),
314                             static_cast<float>(tmp[2]));
315                 glUniform3i(shaders.bond_shader.pbc_cell, x, y, z);
316                 glDrawArraysInstanced(GL_TRIANGLES, 0,
317                                       bond_model_npoly,
318                                       static_cast<GLsizei>(bond_buffer.size()));
319             }
320         }
321     }
322     // cell
323     if(drawCell && curStep->hasCell()){
324         glBindVertexArray(objects.cell_vao);
325         glUseProgram(shaders.cell_shader.program);
326         for(int x=0;x<mult[0];++x){
327             for(int y=0;y<mult[1];++y){
328                 for(int z=0;z<mult[2];++z){
329                     tmp = (off + x*cv[0] + y*cv[1] + z*cv[2]);
330                     glUniform3f(shaders.cell_shader.offset,
331                                 static_cast<float>(tmp[0]),
332                                 static_cast<float>(tmp[1]),
333                                 static_cast<float>(tmp[2]));
334                     glDrawElements(GL_LINES, 24, GL_UNSIGNED_SHORT, nullptr);
335                 }
336             }
337         }
338     }
339 }
340 
drawSel(Vec off,const PBCVec & mult,void * context)341 void GUI::StepData::drawSel(Vec off, const PBCVec &mult, void *context)
342 {
343     // draw Atoms (as setup in atom_vao, shader locations must match)
344     // with selection shader -> color by gl_InstanceID
345     // to seperate Framebuffer
346     const auto &sel_shader = shader_map[context].sel_shader;
347     glEnable(GL_CULL_FACE);
348     glEnable(GL_DEPTH_TEST);
349     glClearColor(1,1,0,0);
350     glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
351     glBindVertexArray(object_map[context].sel_vao);
352     glUseProgram(sel_shader.program);
353     glUniform1f(sel_shader.scale_fac, atRadFac);
354     glUniformMatrix3fv(sel_shader.pos_scale, 1, 0, cell_mat.data());
355     if(curStep->hasCell()){
356         Vec tmp;
357         Mat cv = curStep->getCellVec() * curStep->getCellDim(AtomFmt::Bohr);
358         off -= (mult[0]-1)*cv[0]/2.;
359         off -= (mult[1]-1)*cv[1]/2.;
360         off -= (mult[2]-1)*cv[2]/2.;
361         for(unsigned int x=0;x<mult[0];++x){
362             for(unsigned int y=0;y<mult[1];++y){
363                 for(unsigned int z=0;z<mult[2];++z){
364                     tmp = (off + x*cv[0] + y*cv[1] + z*cv[2]);
365                     glUniform3f(sel_shader.offset,
366                                 static_cast<float>(tmp[0]),
367                                 static_cast<float>(tmp[1]),
368                                 static_cast<float>(tmp[2]));
369                     glUniform1ui(sel_shader.pbc_instance,
370                                  1 + x + y*mult[0] + z*mult[0]*mult[1]);
371                     glDrawArraysInstanced(GL_TRIANGLES, 0,
372                                           atom_model_npoly,
373                                           static_cast<GLsizei>(atom_buffer.size()));
374                 }
375             }
376         }
377     }else{
378         glUniform3f(sel_shader.offset,
379                     static_cast<float>(off[0]),
380                     static_cast<float>(off[1]),
381                     static_cast<float>(off[2]));
382         glUniform1ui(sel_shader.pbc_instance, 1);
383         glDrawArraysInstanced(GL_TRIANGLES, 0,
384                               atom_model_npoly,
385                               static_cast<GLsizei>(atom_buffer.size()));
386     }
387 }
388 
updateGL(void * context)389 void GUI::StepData::updateGL(void *context)
390 {
391     //TODO: separate data-handling somehow
392     auto &objects = object_map[context];
393     // ATOMS
394     auto nat = atom_buffer.size();
395     glBindBuffer(GL_ARRAY_BUFFER, objects.atom_vbo);
396     if (nat != 0u) {
397         glBufferData(GL_ARRAY_BUFFER, static_cast<GLsizeiptr>(nat*sizeof(AtomProp)),
398                      static_cast<void*>(atom_buffer.data()), GL_STREAM_DRAW);
399     } else {
400         glBufferData(GL_ARRAY_BUFFER, 0, nullptr, GL_STREAM_DRAW);
401     }
402     // BONDS
403     glBindBuffer(GL_ARRAY_BUFFER, objects.bond_vbo);
404     if (!bond_buffer.empty()) {
405         glBufferData(GL_ARRAY_BUFFER, static_cast<GLsizeiptr>(bond_buffer.size()*sizeof(BondProp)),
406                      static_cast<void*>(bond_buffer.data()), GL_STREAM_DRAW);
407     } else {
408         glBufferData(GL_ARRAY_BUFFER, 0, nullptr, GL_STREAM_DRAW);
409     }
410     // CELL
411     glBindBuffer(GL_ARRAY_BUFFER, objects.cell_vbo);
412     glBufferData(GL_ARRAY_BUFFER, sizeof(cell_buffer),
413                  static_cast<void*>(cell_buffer.data()), GL_STREAM_DRAW);
414 }
415 
416 // TODO: radii should probably be factored in in shader
update(Step * step,bool useVdW,float atRadFac,float bondRad)417 void GUI::StepData::update(Step* step,
418                            bool useVdW, float atRadFac,
419                            float bondRad)
420 {
421     curStep = step;
422     for(auto &[context, state]: instance_map){
423         state.synchronized = false;
424     }
425     this->atRadFac = atRadFac;
426 
427 // CELL
428     Mat cv = curStep->getCellVec() * curStep->getCellDim(AtomFmt::Bohr);
429     cell_buffer = {0,0,0,
430                   static_cast<float>(cv[0][0]), static_cast<float>(cv[0][1]), static_cast<float>(cv[0][2]),
431                   static_cast<float>(cv[1][0]), static_cast<float>(cv[1][1]), static_cast<float>(cv[1][2]),
432                   static_cast<float>(cv[2][0]), static_cast<float>(cv[2][1]), static_cast<float>(cv[2][2]),
433                   static_cast<float>(cv[0][0]+cv[1][0]), static_cast<float>(cv[0][1]+cv[1][1]), static_cast<float>(cv[0][2]+cv[1][2]),
434                   static_cast<float>(cv[0][0]+cv[2][0]), static_cast<float>(cv[0][1]+cv[2][1]), static_cast<float>(cv[0][2]+cv[2][2]),
435                   static_cast<float>(cv[1][0]+cv[2][0]), static_cast<float>(cv[1][1]+cv[2][1]), static_cast<float>(cv[1][2]+cv[2][2]),
436                   static_cast<float>(cv[0][0]+cv[1][0]+cv[2][0]), static_cast<float>(cv[0][1]+cv[1][1]+cv[2][1]), static_cast<float>(cv[0][2]+cv[1][2]+cv[2][2]),
437                   };
438     Mat tmp_mat;
439     if(curStep->getFmt() == AtomFmt::Crystal){
440         tmp_mat = curStep->getCellVec();
441     }else{
442         tmp_mat = {{{{1,0,0}}, {{0,1,0}}, {{0,0,1}}}};
443     }
444     switch(curStep->getFmt()){
445     case AtomFmt::Angstrom:
446         tmp_mat *= Vipster::invbohr;
447         break;
448     case AtomFmt::Crystal:
449     case AtomFmt::Alat:
450         tmp_mat *= curStep->getCellDim(AtomFmt::Bohr);
451         break;
452     default:
453         break;
454     }
455     cell_mat = {{static_cast<float>(tmp_mat[0][0]),
456                  static_cast<float>(tmp_mat[1][0]),
457                  static_cast<float>(tmp_mat[2][0]),
458                  static_cast<float>(tmp_mat[0][1]),
459                  static_cast<float>(tmp_mat[1][1]),
460                  static_cast<float>(tmp_mat[2][1]),
461                  static_cast<float>(tmp_mat[0][2]),
462                  static_cast<float>(tmp_mat[1][2]),
463                  static_cast<float>(tmp_mat[2][2])}};
464 
465 // ATOMS
466     atom_buffer.clear();
467     atom_buffer.reserve(curStep->getNat());
468     if(useVdW){
469         for (const auto& at: *curStep){
470             atom_buffer.push_back({{static_cast<float>(at.coord[0]),
471                                     static_cast<float>(at.coord[1]),
472                                     static_cast<float>(at.coord[2])},
473                                    static_cast<float>(at.type->vdwr), at.type->col,
474                                    static_cast<uint8_t>(at.properties->flags[AtomProperties::Hidden])});
475         }
476     }else{
477         for (const auto& at: *curStep){
478             atom_buffer.push_back({{static_cast<float>(at.coord[0]),
479                                     static_cast<float>(at.coord[1]),
480                                     static_cast<float>(at.coord[2])},
481                                    static_cast<float>(at.type->covr), at.type->col,
482                                    static_cast<uint8_t>(at.properties->flags[AtomProperties::Hidden])});
483         }
484     }
485 
486 // BONDS
487     constexpr Vec x_axis{{1,0,0}};
488     const auto& bonds = curStep->getBonds();
489     const auto& elements = curStep->getAtoms().elements;
490     const auto& at_coord = curStep->getAtoms().coordinates;
491     const auto& at_prop = curStep->getAtoms().properties;
492     auto fmt = curStep->getFmt();
493     auto fmt_ax_crystal = cv;
494     float c, s, ic;
495     bond_buffer.clear();
496     bond_buffer.reserve(bonds.size());
497     switch(fmt){
498     case AtomFmt::Crystal:
499         cv = Mat{{{{1,0,0}}, {{0,1,0}}, {{0,0,1}}}};
500         break;
501     case AtomFmt::Alat:
502         cv = curStep->getCellVec();
503         break;
504     case AtomFmt::Angstrom:
505         cv = curStep->getCellVec() * curStep->getCellDim(AtomFmt::Angstrom);
506         break;
507     case AtomFmt::Bohr:
508         break;
509     }
510     for(const Bond& bd:bonds){
511         auto at_pos1 = at_coord[bd.at1];
512         auto at_pos2 = at_coord[bd.at2];
513         if (bd.diff[0]>0)     { at_pos2 += bd.diff[0]*cv[0]; }
514         else if (bd.diff[0]<0){ at_pos1 -= bd.diff[0]*cv[0]; }
515         if (bd.diff[1]>0)     { at_pos2 += bd.diff[1]*cv[1]; }
516         else if (bd.diff[1]<0){ at_pos1 -= bd.diff[1]*cv[1]; }
517         if (bd.diff[2]>0)     { at_pos2 += bd.diff[2]*cv[2]; }
518         else if (bd.diff[2]<0){ at_pos1 -= bd.diff[2]*cv[2]; }
519         auto bond_axis = at_pos1 - at_pos2;
520         if(fmt == AtomFmt::Crystal){
521             bond_axis = bond_axis * fmt_ax_crystal;
522         }
523         auto bond_pos = (at_pos1+at_pos2)/2;
524         const auto& col1 = bd.type ? bd.type->second : elements[bd.at1]->second.col;
525         const auto& col2 = bd.type ? bd.type->second : elements[bd.at2]->second.col;
526         // handle bonds parallel to x-axis
527         if(std::abs(bond_axis[1])<std::numeric_limits<Vec::value_type>::epsilon()&&
528            std::abs(bond_axis[2])<std::numeric_limits<Vec::value_type>::epsilon()){
529             c = std::copysign(1.f, static_cast<float>(bond_axis[0]));
530             bond_buffer.push_back({
531                 //mat3 with rotation and scaling
532                 {static_cast<float>(bd.dist)*c, 0., 0.,
533                  0., bondRad, 0.,
534                  0., 0., bondRad*c},
535                 //vec3 with position in modelspace
536                 {static_cast<float>(bond_pos[0]),
537                  static_cast<float>(bond_pos[1]),
538                  static_cast<float>(bond_pos[2])},
539                 //faux uvec4 with render criteria
540                 {static_cast<int16_t>(std::abs(bd.diff[0])),
541                  static_cast<int16_t>(std::abs(bd.diff[1])),
542                  static_cast<int16_t>(std::abs(bd.diff[2])),
543                  static_cast<int16_t>(at_prop[bd.at1].flags[AtomProperties::Hidden] ||
544                                       at_prop[bd.at2].flags[AtomProperties::Hidden])},
545                  col1, col2});
546         }else{
547             // all other bonds
548             auto rot_axis = -Vec_cross(bond_axis, x_axis);
549             rot_axis /= Vec_length(rot_axis);
550             c = static_cast<float>(Vec_dot(bond_axis, x_axis)/Vec_length(bond_axis));
551             ic = 1-c;
552             s = -std::sqrt(1-c*c);
553             const auto dist = static_cast<float>(bd.dist);
554             const float ax[3] = {static_cast<float>(rot_axis[0]),
555                                  static_cast<float>(rot_axis[1]),
556                                  static_cast<float>(rot_axis[2])};
557             bond_buffer.push_back({
558                 //mat3 with rotation and scaling
559                 {dist*(ic*ax[0]*ax[0]+c),
560                  dist*(ic*ax[0]*ax[1]-s*ax[2]),
561                  dist*(ic*ax[0]*ax[2]+s*ax[1]),
562                  bondRad*(ic*ax[1]*ax[0]+s*ax[2]),
563                  bondRad*(ic*ax[1]*ax[1]+c),
564                  bondRad*(ic*ax[1]*ax[2]-s*ax[0]),
565                  bondRad*(ic*ax[2]*ax[0]-s*ax[1]),
566                  bondRad*(ic*ax[2]*ax[1]+s*ax[0]),
567                  bondRad*(ic*ax[2]*ax[2]+c)},
568                 //vec3 with position in modelspace
569                 {static_cast<float>(bond_pos[0]),
570                  static_cast<float>(bond_pos[1]),
571                  static_cast<float>(bond_pos[2])},
572                 //faux uvec4 with render criteria
573                 {static_cast<int16_t>(std::abs(bd.diff[0])),
574                  static_cast<int16_t>(std::abs(bd.diff[1])),
575                  static_cast<int16_t>(std::abs(bd.diff[2])),
576                  static_cast<int16_t>(at_prop[bd.at1].flags[AtomProperties::Hidden] ||
577                                       at_prop[bd.at2].flags[AtomProperties::Hidden])},
578                 //2*vec4 with colors
579                 col1, col2});
580         }
581     }
582 }
583