1 /* molecule, Copyright (c) 2001-2012 Jamie Zawinski <jwz@jwz.org>
2 * Draws molecules, based on coordinates from PDB (Protein Data Base) files.
4 * Permission to use, copy, modify, distribute, and sell this software and its
5 * documentation for any purpose is hereby granted without fee, provided that
6 * the above copyright notice appear in all copies and that both that
7 * copyright notice and this permission notice appear in supporting
8 * documentation. No representations are made about the suitability of this
9 * software for any purpose. It is provided "as is" without express or
14 /* Documentation on the PDB file format:
15 http://en.wikipedia.org/wiki/Protein_Data_Bank_%28file_format%29
16 http://www.wwpdb.org/docs.html
17 http://www.wwpdb.org/documentation/format32/v3.2.html
18 http://www.wwpdb.org/documentation/format32/sect9.html
19 http://www.rcsb.org/pdb/file_formats/pdb/pdbguide2.2/guide2.2_frame.html
21 Good source of PDB files:
22 http://www.sci.ouc.bc.ca/chem/molecule/molecule.html
23 http://www.umass.edu/microbio/rasmol/whereget.htm
24 http://www.wwpdb.org/docs.html
28 # define ATOM_FONT "-*-helvetica-medium-r-normal-*-180-*"
30 # define ATOM_FONT "-*-helvetica-medium-r-normal-*-240-*"
33 #define DEFAULTS "*delay: 10000 \n" \
34 "*showFPS: False \n" \
35 "*wireframe: False \n" \
36 "*atomFont: " ATOM_FONT "\n" \
37 "*atomFont2: -*-helvetica-bold-r-normal-*-80-*\n" \
38 "*titleFont: -*-helvetica-medium-r-normal-*-180-*\n" \
39 "*noLabelThreshold: 30 \n" \
40 "*wireframeThreshold: 150 \n" \
42 # define refresh_molecule 0
43 # define release_molecule 0
45 #define countof(x) (sizeof((x))/sizeof((*x)))
47 #include "xlockmore.h"
53 #include "gltrackball.h"
55 #ifdef USE_GL /* whole file */
57 #include <sys/types.h>
62 #define DEF_TIMEOUT "20"
63 #define DEF_SPIN "XYZ"
64 #define DEF_WANDER "False"
65 #define DEF_LABELS "True"
66 #define DEF_TITLES "True"
67 #define DEF_ATOMS "True"
68 #define DEF_BONDS "True"
69 #define DEF_ESHELLS "True"
70 #define DEF_BBOX "False"
71 #define DEF_SHELL_ALPHA "0.3"
72 #define DEF_MOLECULE "(default)"
73 #define DEF_VERBOSE "False"
75 #define SPHERE_SLICES 48 /* how densely to render spheres */
76 #define SPHERE_STACKS 24
78 #define SMOOTH_TUBE /* whether to have smooth or faceted tubes */
81 # define TUBE_FACES 12 /* how densely to render tubes */
86 #define SPHERE_SLICES_2 14
87 #define SPHERE_STACKS_2 8
88 #define TUBE_FACES_2 6
92 __extension__ /* don't warn about "string length is greater than the length
93 ISO C89 compilers are required to support" when includng
94 the following data file... */
96 static const char * const builtin_pdb_data[] = {
97 # include "molecules.h"
110 const char *text_color;
115 /* These are the traditional colors used to render these atoms,
116 and their approximate size in angstroms.
118 static const atom_data all_atom_data[] = {
119 { "H", 1.17, 0.40, "#FFFFFF", "#000000", { 0, }},
120 { "C", 1.75, 0.58, "#999999", "#FFFFFF", { 0, }},
121 { "CA", 1.80, 0.60, "#0000FF", "#ADD8E6", { 0, }},
122 { "N", 1.55, 0.52, "#A2B5CD", "#EE99FF", { 0, }},
123 { "O", 1.40, 0.47, "#FF0000", "#FFB6C1", { 0, }},
124 { "P", 1.28, 0.43, "#9370DB", "#DB7093", { 0, }},
125 { "S", 1.80, 0.60, "#8B8B00", "#FFFF00", { 0, }},
126 { "bond", 0, 0, "#B3B3B3", "#FFFF00", { 0, }},
127 { "*", 1.40, 0.47, "#008B00", "#90EE90", { 0, }}
132 int id; /* sequence number in the PDB file */
133 const char *label; /* The atom name */
134 GLfloat x, y, z; /* position in 3-space (angstroms) */
135 const atom_data *data; /* computed: which style of atom this is */
139 int from, to; /* atom sequence numbers */
140 int strength; /* how many bonds are between these two atoms */
145 const char *label; /* description of this compound */
146 int natoms, atoms_size;
147 int nbonds, bonds_size;
148 molecule_atom *atoms;
149 molecule_bond *bonds;
154 GLXContext *glx_context;
156 trackball_state *trackball;
159 GLfloat molecule_size; /* max dimension of molecule bounding box */
161 GLfloat no_label_threshold; /* Things happen when molecules are huge */
162 GLfloat wireframe_threshold;
164 int which; /* which of the molecules is being shown */
168 int mode; /* 0 = normal, 1 = out, 2 = in */
171 GLuint molecule_dlist;
174 # ifdef HAVE_GLBITMAP
175 XFontStruct *xfont1, *xfont2, *xfont3;
176 GLuint font1_dlist, font2_dlist, font3_dlist;
178 texture_font_data *font1_data, *font2_data, *font3_data;
189 } molecule_configuration;
192 static molecule_configuration *mcs = NULL;
195 static char *molecule_str;
196 static char *do_spin;
197 static Bool do_wander;
198 static Bool do_titles;
199 static Bool do_labels;
200 static Bool do_atoms;
201 static Bool do_bonds;
202 static Bool do_shells;
204 static Bool verbose_p;
205 static GLfloat shell_alpha;
208 static Bool orig_do_labels, orig_do_atoms, orig_do_bonds, orig_do_shells,
212 static XrmOptionDescRec opts[] = {
213 { "-molecule", ".molecule", XrmoptionSepArg, 0 },
214 { "-timeout", ".timeout", XrmoptionSepArg, 0 },
215 { "-spin", ".spin", XrmoptionSepArg, 0 },
216 { "+spin", ".spin", XrmoptionNoArg, "" },
217 { "-wander", ".wander", XrmoptionNoArg, "True" },
218 { "+wander", ".wander", XrmoptionNoArg, "False" },
219 { "-labels", ".labels", XrmoptionNoArg, "True" },
220 { "+labels", ".labels", XrmoptionNoArg, "False" },
221 { "-titles", ".titles", XrmoptionNoArg, "True" },
222 { "+titles", ".titles", XrmoptionNoArg, "False" },
223 { "-atoms", ".atoms", XrmoptionNoArg, "True" },
224 { "+atoms", ".atoms", XrmoptionNoArg, "False" },
225 { "-bonds", ".bonds", XrmoptionNoArg, "True" },
226 { "+bonds", ".bonds", XrmoptionNoArg, "False" },
227 { "-shells", ".eshells", XrmoptionNoArg, "True" },
228 { "+shells", ".eshells", XrmoptionNoArg, "False" },
229 { "-shell-alpha", ".shellAlpha", XrmoptionSepArg, 0 },
230 { "-bbox", ".bbox", XrmoptionNoArg, "True" },
231 { "+bbox", ".bbox", XrmoptionNoArg, "False" },
232 { "-verbose", ".verbose", XrmoptionNoArg, "True" },
235 static argtype vars[] = {
236 {&molecule_str, "molecule", "Molecule", DEF_MOLECULE, t_String},
237 {&timeout, "timeout", "Seconds", DEF_TIMEOUT, t_Int},
238 {&do_spin, "spin", "Spin", DEF_SPIN, t_String},
239 {&do_wander, "wander", "Wander", DEF_WANDER, t_Bool},
240 {&do_atoms, "atoms", "Atoms", DEF_ATOMS, t_Bool},
241 {&do_bonds, "bonds", "Bonds", DEF_BONDS, t_Bool},
242 {&do_shells, "eshells", "EShells", DEF_ESHELLS, t_Bool},
243 {&do_labels, "labels", "Labels", DEF_LABELS, t_Bool},
244 {&do_titles, "titles", "Titles", DEF_TITLES, t_Bool},
245 {&do_bbox, "bbox", "BBox", DEF_BBOX, t_Bool},
246 {&shell_alpha, "shellAlpha", "ShellAlpha", DEF_SHELL_ALPHA, t_Float},
247 {&verbose_p, "verbose", "Verbose", DEF_VERBOSE, t_Bool},
250 ENTRYPOINT ModeSpecOpt molecule_opts = {countof(opts), opts, countof(vars), vars, NULL};
258 sphere (molecule_configuration *mc,
259 GLfloat x, GLfloat y, GLfloat z, GLfloat diameter, Bool wire)
261 int stacks = (mc->scale_down ? SPHERE_STACKS_2 : SPHERE_STACKS);
262 int slices = (mc->scale_down ? SPHERE_SLICES_2 : SPHERE_SLICES);
265 glTranslatef (x, y, z);
266 glScalef (diameter, diameter, diameter);
267 unit_sphere (stacks, slices, wire);
270 return stacks * slices;
275 load_fonts (ModeInfo *mi)
277 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
278 # ifdef HAVE_GLBITMAP
279 load_font (mi->dpy, "atomFont", &mc->xfont1, &mc->font1_dlist);
280 load_font (mi->dpy, "atomFont2", &mc->xfont2, &mc->font2_dlist);
281 load_font (mi->dpy, "titleFont", &mc->xfont3, &mc->font3_dlist);
283 mc->font1_data = load_texture_font (mi->dpy, "atomFont");
284 mc->font2_data = load_texture_font (mi->dpy, "atomFont2");
285 mc->font3_data = load_texture_font (mi->dpy, "titleFont");
290 static const atom_data *
291 get_atom_data (const char *atom_name)
294 const atom_data *d = 0;
295 char *n = strdup (atom_name);
299 while (!isalpha(*n)) n++;
301 while (L > 0 && !isalpha(n[L-1]))
304 for (i = 0; i < countof(all_atom_data); i++)
306 d = &all_atom_data[i];
307 if (!strcasecmp (n, all_atom_data[i].name))
317 set_atom_color (ModeInfo *mi, const molecule_atom *a,
318 Bool font_p, GLfloat alpha)
326 d = get_atom_data ("bond");
330 gl_color[0] = d->gl_color[4];
331 gl_color[1] = d->gl_color[5];
332 gl_color[2] = d->gl_color[6];
333 gl_color[3] = d->gl_color[7];
337 gl_color[0] = d->gl_color[0];
338 gl_color[1] = d->gl_color[1];
339 gl_color[2] = d->gl_color[2];
340 gl_color[3] = d->gl_color[3];
343 if (gl_color[3] == 0)
345 const char *string = !font_p ? d->color : d->text_color;
347 if (!XParseColor (mi->dpy, mi->xgwa.colormap, string, &xcolor))
349 fprintf (stderr, "%s: unparsable color in %s: %s\n", progname,
350 (a ? a->label : d->name), string);
354 gl_color[0] = xcolor.red / 65536.0;
355 gl_color[1] = xcolor.green / 65536.0;
356 gl_color[2] = xcolor.blue / 65536.0;
361 /* If we're not drawing atoms, and the color is black, use white instead.
362 This is a kludge so that H can have black text over its white ball,
363 but the text still shows up if balls are off.
365 if (font_p && !do_atoms &&
366 gl_color[0] == 0 && gl_color[1] == 0 && gl_color[2] == 0)
368 gl_color[0] = gl_color[1] = gl_color[2] = 1;
372 glColor4f (gl_color[0], gl_color[1], gl_color[2], gl_color[3]);
374 glMaterialfv (GL_FRONT, GL_AMBIENT_AND_DIFFUSE, gl_color);
379 atom_size (const molecule_atom *a)
382 return a->data->size2;
384 return a->data->size;
388 static molecule_atom *
389 get_atom (molecule_atom *atoms, int natoms, int id)
393 /* quick short-circuit */
396 if (atoms[id].id == id)
398 if (id > 0 && atoms[id-1].id == id)
400 if (id < natoms-1 && atoms[id+1].id == id)
404 for (i = 0; i < natoms; i++)
405 if (id == atoms[i].id)
408 fprintf (stderr, "%s: no atom %d\n", progname, id);
414 molecule_bounding_box (ModeInfo *mi,
415 GLfloat *x1, GLfloat *y1, GLfloat *z1,
416 GLfloat *x2, GLfloat *y2, GLfloat *z2)
418 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
419 molecule *m = &mc->molecules[mc->which];
424 *x1 = *y1 = *z1 = *x2 = *y2 = *z2 = 0;
428 *x1 = *x2 = m->atoms[0].x;
429 *y1 = *y2 = m->atoms[0].y;
430 *z1 = *z2 = m->atoms[0].z;
433 for (i = 1; i < m->natoms; i++)
435 if (m->atoms[i].x < *x1) *x1 = m->atoms[i].x;
436 if (m->atoms[i].y < *y1) *y1 = m->atoms[i].y;
437 if (m->atoms[i].z < *z1) *z1 = m->atoms[i].z;
439 if (m->atoms[i].x > *x2) *x2 = m->atoms[i].x;
440 if (m->atoms[i].y > *y2) *y2 = m->atoms[i].y;
441 if (m->atoms[i].z > *z2) *z2 = m->atoms[i].z;
454 draw_bounding_box (ModeInfo *mi)
456 static const GLfloat c1[4] = { 0.2, 0.2, 0.4, 1.0 };
457 static const GLfloat c2[4] = { 1.0, 0.0, 0.0, 1.0 };
458 int wire = MI_IS_WIREFRAME(mi);
459 GLfloat x1, y1, z1, x2, y2, z2;
460 molecule_bounding_box (mi, &x1, &y1, &z1, &x2, &y2, &z2);
462 glMaterialfv (GL_FRONT, GL_AMBIENT_AND_DIFFUSE, c1);
465 glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
467 glVertex3f(x1, y1, z1); glVertex3f(x1, y1, z2);
468 glVertex3f(x2, y1, z2); glVertex3f(x2, y1, z1);
470 glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
471 glNormal3f(0, -1, 0);
472 glVertex3f(x2, y2, z1); glVertex3f(x2, y2, z2);
473 glVertex3f(x1, y2, z2); glVertex3f(x1, y2, z1);
475 glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
477 glVertex3f(x1, y1, z1); glVertex3f(x2, y1, z1);
478 glVertex3f(x2, y2, z1); glVertex3f(x1, y2, z1);
480 glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
481 glNormal3f(0, 0, -1);
482 glVertex3f(x1, y2, z2); glVertex3f(x2, y2, z2);
483 glVertex3f(x2, y1, z2); glVertex3f(x1, y1, z2);
485 glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
487 glVertex3f(x1, y2, z1); glVertex3f(x1, y2, z2);
488 glVertex3f(x1, y1, z2); glVertex3f(x1, y1, z1);
490 glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
491 glNormal3f(-1, 0, 0);
492 glVertex3f(x2, y1, z1); glVertex3f(x2, y1, z2);
493 glVertex3f(x2, y2, z2); glVertex3f(x2, y2, z1);
496 glDisable (GL_LIGHTING);
498 glColor3f (c2[0], c2[1], c2[2]);
500 if (x1 > 0) x1 = 0; if (x2 < 0) x2 = 0;
501 if (y1 > 0) y1 = 0; if (y2 < 0) y2 = 0;
502 if (z1 > 0) z1 = 0; if (z2 < 0) z2 = 0;
503 glVertex3f(x1, 0, 0); glVertex3f(x2, 0, 0);
504 glVertex3f(0 , y1, 0); glVertex3f(0, y2, 0);
505 glVertex3f(0, 0, z1); glVertex3f(0, 0, z2);
509 glEnable (GL_LIGHTING);
513 /* Since PDB files don't always have the molecule centered around the
514 origin, and since some molecules are pretty large, scale and/or
515 translate so that the whole molecule is visible in the window.
518 ensure_bounding_box_visible (ModeInfo *mi)
520 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
522 GLfloat x1, y1, z1, x2, y2, z2;
525 GLfloat max_size = 10; /* don't bother scaling down if the molecule
526 is already smaller than this */
528 molecule_bounding_box (mi, &x1, &y1, &z1, &x2, &y2, &z2);
533 size = (w > h ? w : h);
534 size = (size > d ? size : d);
536 mc->molecule_size = size;
542 GLfloat scale = max_size / size;
543 glScalef (scale, scale, scale);
545 mc->scale_down = scale < 0.3;
548 glTranslatef (-(x1 + w/2),
555 /* Constructs the GL shapes of the current molecule
558 build_molecule (ModeInfo *mi, Bool transparent_p)
560 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
561 int wire = MI_IS_WIREFRAME(mi);
563 GLfloat alpha = transparent_p ? shell_alpha : 1.0;
566 molecule *m = &mc->molecules[mc->which];
570 glDisable(GL_CULL_FACE);
571 glDisable(GL_LIGHTING);
572 glDisable(GL_LIGHT0);
573 glDisable(GL_DEPTH_TEST);
574 glDisable(GL_NORMALIZE);
575 glDisable(GL_CULL_FACE);
579 glEnable(GL_CULL_FACE);
580 glEnable(GL_LIGHTING);
582 glEnable(GL_DEPTH_TEST);
583 glEnable(GL_NORMALIZE);
584 glEnable(GL_CULL_FACE);
588 set_atom_color (mi, 0, False, alpha);
591 for (i = 0; i < m->nbonds; i++)
593 const molecule_bond *b = &m->bonds[i];
594 const molecule_atom *from = get_atom (m->atoms, m->natoms, b->from);
595 const molecule_atom *to = get_atom (m->atoms, m->natoms, b->to);
600 glVertex3f(from->x, from->y, from->z);
601 glVertex3f(to->x, to->y, to->z);
607 int faces = (mc->scale_down ? TUBE_FACES_2 : TUBE_FACES);
613 GLfloat thickness = 0.07 * b->strength;
614 GLfloat cap_size = 0.03;
618 polys += tube (from->x, from->y, from->z,
621 faces, smooth, (!do_atoms || do_shells), wire);
625 if (!wire && do_atoms)
626 for (i = 0; i < m->natoms; i++)
628 const molecule_atom *a = &m->atoms[i];
629 GLfloat size = atom_size (a);
630 set_atom_color (mi, a, False, alpha);
631 polys += sphere (mc, a->x, a->y, a->z, size, wire);
634 if (do_bbox && !transparent_p)
636 draw_bounding_box (mi);
640 mc->polygon_count += polys;
648 push_atom (molecule *m,
649 int id, const char *label,
650 GLfloat x, GLfloat y, GLfloat z)
653 if (m->atoms_size < m->natoms)
656 m->atoms = (molecule_atom *) realloc (m->atoms,
657 m->atoms_size * sizeof(*m->atoms));
659 m->atoms[m->natoms-1].id = id;
660 m->atoms[m->natoms-1].label = label;
661 m->atoms[m->natoms-1].x = x;
662 m->atoms[m->natoms-1].y = y;
663 m->atoms[m->natoms-1].z = z;
664 m->atoms[m->natoms-1].data = get_atom_data (label);
669 push_bond (molecule *m, int from, int to)
673 for (i = 0; i < m->nbonds; i++)
674 if ((m->bonds[i].from == from && m->bonds[i].to == to) ||
675 (m->bonds[i].to == from && m->bonds[i].from == to))
677 m->bonds[i].strength++;
682 if (m->bonds_size < m->nbonds)
685 m->bonds = (molecule_bond *) realloc (m->bonds,
686 m->bonds_size * sizeof(*m->bonds));
688 m->bonds[m->nbonds-1].from = from;
689 m->bonds[m->nbonds-1].to = to;
690 m->bonds[m->nbonds-1].strength = 1;
695 parse_error (const char *file, int lineno, const char *line)
697 fprintf (stderr, "%s: %s: parse error, line %d: %s\n",
698 progname, file, lineno, line);
703 /* This function is crap.
706 parse_pdb_data (molecule *m, const char *data, const char *filename, int line)
708 const char *s = data;
712 if ((!m->label || !*m->label) &&
713 (!strncmp (s, "HEADER", 6) || !strncmp (s, "COMPND", 6)))
715 char *name = calloc (1, 100);
722 while (isspace(*n2)) n2++;
724 ss = strchr (n2, '\n');
726 ss = strchr (n2, '\r');
729 ss = n2+strlen(n2)-1;
730 while (isspace(*ss) && ss > n2)
733 if (strlen (n2) > 4 &&
734 !strcmp (n2 + strlen(n2) - 4, ".pdb"))
735 n2[strlen(n2)-4] = 0;
737 if (m->label) free ((char *) m->label);
738 m->label = strdup (n2);
741 else if (!strncmp (s, "TITLE ", 6) ||
742 !strncmp (s, "HEADER", 6) ||
743 !strncmp (s, "COMPND", 6) ||
744 !strncmp (s, "AUTHOR", 6) ||
745 !strncmp (s, "REVDAT", 6) ||
746 !strncmp (s, "SOURCE", 6) ||
747 !strncmp (s, "EXPDTA", 6) ||
748 !strncmp (s, "JRNL ", 6) ||
749 !strncmp (s, "REMARK", 6) ||
750 !strncmp (s, "SEQRES", 6) ||
751 !strncmp (s, "HET ", 6) ||
752 !strncmp (s, "FORMUL", 6) ||
753 !strncmp (s, "CRYST1", 6) ||
754 !strncmp (s, "ORIGX1", 6) ||
755 !strncmp (s, "ORIGX2", 6) ||
756 !strncmp (s, "ORIGX3", 6) ||
757 !strncmp (s, "SCALE1", 6) ||
758 !strncmp (s, "SCALE2", 6) ||
759 !strncmp (s, "SCALE3", 6) ||
760 !strncmp (s, "MASTER", 6) ||
761 !strncmp (s, "KEYWDS", 6) ||
762 !strncmp (s, "DBREF ", 6) ||
763 !strncmp (s, "HETNAM", 6) ||
764 !strncmp (s, "HETSYN", 6) ||
765 !strncmp (s, "HELIX ", 6) ||
766 !strncmp (s, "LINK ", 6) ||
767 !strncmp (s, "MTRIX1", 6) ||
768 !strncmp (s, "MTRIX2", 6) ||
769 !strncmp (s, "MTRIX3", 6) ||
770 !strncmp (s, "SHEET ", 6) ||
771 !strncmp (s, "CISPEP", 6) ||
773 !strncmp (s, "SEQADV", 6) ||
774 !strncmp (s, "SITE ", 5) ||
775 !strncmp (s, "FTNOTE", 6) ||
776 !strncmp (s, "MODEL ", 5) ||
777 !strncmp (s, "ENDMDL", 6) ||
778 !strncmp (s, "SPRSDE", 6) ||
779 !strncmp (s, "MODRES", 6) ||
781 !strncmp (s, "GENERATED BY", 12) ||
782 !strncmp (s, "TER ", 4) ||
783 !strncmp (s, "END ", 4) ||
784 !strncmp (s, "TER\n", 4) ||
785 !strncmp (s, "END\n", 4) ||
786 !strncmp (s, "\n", 1))
789 else if (!strncmp (s, "ATOM ", 7))
792 const char *end = strchr (s, '\n');
794 char *name = (char *) calloc (1, 4);
795 GLfloat x = -999, y = -999, z = -999;
797 if (1 != sscanf (s+7, " %d ", &id))
798 parse_error (filename, line, s);
800 /* Use the "atom name" field if that is all that is available. */
801 strncpy (name, s+12, 3);
803 /* But prefer the "element" field. */
804 if (L > 77 && !isspace(s[77])) {
805 /* fprintf(stderr, " \"%s\" -> ", name); */
809 /* fprintf(stderr, "\"%s\"\n", name); */
812 while (isspace(*name)) name++;
813 ss = name + strlen(name)-1;
814 while (isspace(*ss) && ss > name)
822 if (3 != sscanf (s + 32, " %f %f %f ", &x, &y, &z))
823 parse_error (filename, line, s);
826 fprintf (stderr, "%s: %s: %d: atom: %d \"%s\" %9.4f %9.4f %9.4f\n",
827 progname, filename, line,
830 push_atom (m, id, name, x, y, z);
832 else if (!strncmp (s, "HETATM ", 7))
835 char *name = (char *) calloc (1, 4);
836 GLfloat x = -999, y = -999, z = -999;
838 if (1 != sscanf (s+7, " %d ", &id))
839 parse_error (filename, line, s);
841 strncpy (name, s+12, 3);
842 while (isspace(*name)) name++;
843 ss = name + strlen(name)-1;
844 while (isspace(*ss) && ss > name)
846 if (3 != sscanf (s + 30, " %f %f %f ", &x, &y, &z))
847 parse_error (filename, line, s);
849 fprintf (stderr, "%s: %s: %d: atom: %d \"%s\" %9.4f %9.4f %9.4f\n",
850 progname, filename, line,
853 push_atom (m, id, name, x, y, z);
855 else if (!strncmp (s, "CONECT ", 7))
858 int i = sscanf (s + 8, " %d %d %d %d %d %d %d %d %d %d %d %d ",
859 &atoms[0], &atoms[1], &atoms[2], &atoms[3],
860 &atoms[4], &atoms[5], &atoms[6], &atoms[7],
861 &atoms[8], &atoms[9], &atoms[10], &atoms[11]);
863 for (j = 1; j < i; j++)
867 fprintf (stderr, "%s: %s: %d: bond: %d %d\n",
868 progname, filename, line, atoms[0], atoms[j]);
870 push_bond (m, atoms[0], atoms[j]);
875 char *s1 = strdup (s);
876 for (ss = s1; *ss && *ss != '\n'; ss++)
879 fprintf (stderr, "%s: %s: %d: unrecognised line: %s\n",
880 progname, filename, line, s1);
883 while (*s && *s != '\n')
894 parse_pdb_file (molecule *m, const char *name)
897 int buf_size = 40960;
901 in = fopen(name, "r");
904 char *buf = (char *) malloc(1024 + strlen(name));
905 sprintf(buf, "%s: error reading \"%s\"", progname, name);
910 buf = (char *) malloc (buf_size);
912 while (fgets (buf, buf_size-1, in))
915 for (s = buf; *s; s++)
916 if (*s == '\r') *s = '\n';
917 parse_pdb_data (m, buf, name, line++);
925 fprintf (stderr, "%s: file %s contains no atomic coordinates!\n",
930 if (!m->nbonds && do_bonds)
932 fprintf (stderr, "%s: warning: file %s contains no atomic bond info.\n",
939 #endif /* LOAD_FILES */
942 typedef struct { char *atom; int count; } atom_and_count;
944 /* When listing the components of a molecule, the convention is to put the
945 carbon atoms first, the hydrogen atoms second, and the other atom types
946 sorted alphabetically after that (although for some molecules, the usual
947 order is different: we special-case a few of those.)
950 cmp_atoms (const void *aa, const void *bb)
952 const atom_and_count *a = (atom_and_count *) aa;
953 const atom_and_count *b = (atom_and_count *) bb;
954 if (!a->atom) return 1;
955 if (!b->atom) return -1;
956 if (!strcmp(a->atom, "C")) return -1;
957 if (!strcmp(b->atom, "C")) return 1;
958 if (!strcmp(a->atom, "H")) return -1;
959 if (!strcmp(b->atom, "H")) return 1;
960 return strcmp (a->atom, b->atom);
963 static void special_case_formula (char *f);
966 generate_molecule_formula (molecule *m)
968 char *buf = (char *) malloc (m->natoms * 10);
971 atom_and_count counts[200];
972 memset (counts, 0, sizeof(counts));
974 for (i = 0; i < m->natoms; i++)
977 char *a = (char *) m->atoms[i].label;
979 while (!isalpha(*a)) a++;
981 for (e = a; isalpha(*e); e++);
983 while (counts[j].atom && !!strcmp(a, counts[j].atom))
993 while (counts[i].atom) i++;
994 qsort (counts, i, sizeof(*counts), cmp_atoms);
997 while (counts[i].atom)
999 strcat (s, counts[i].atom);
1000 free (counts[i].atom);
1002 if (counts[i].count > 1)
1003 sprintf (s, "[%d]", counts[i].count); /* use [] to get subscripts */
1008 special_case_formula (buf);
1010 if (!m->label) m->label = strdup("");
1011 s = (char *) malloc (strlen (m->label) + strlen (buf) + 2);
1012 strcpy (s, m->label);
1015 free ((char *) m->label);
1020 /* thanks to Rene Uittenbogaard <ruittenb@wish.nl> */
1022 special_case_formula (char *f)
1024 if (!strcmp(f, "H[2]Be")) strcpy(f, "BeH[2]");
1025 else if (!strcmp(f, "H[3]B")) strcpy(f, "BH[3]");
1026 else if (!strcmp(f, "H[3]N")) strcpy(f, "NH[3]");
1027 else if (!strcmp(f, "CHN")) strcpy(f, "HCN");
1028 else if (!strcmp(f, "CKN")) strcpy(f, "KCN");
1029 else if (!strcmp(f, "H[4]N[2]")) strcpy(f, "N[2]H[4]");
1030 else if (!strcmp(f, "Cl[3]P")) strcpy(f, "PCl[3]");
1031 else if (!strcmp(f, "Cl[5]P")) strcpy(f, "PCl[5]");
1036 insert_vertical_whitespace (char *string)
1040 if ((string[0] == ',' ||
1042 string[0] == ':') &&
1044 string[0] = ' ', string[1] = '\n';
1050 /* Construct the molecule data from either: the builtins; or from
1051 the (one) .pdb file specified with -molecule.
1054 load_molecules (ModeInfo *mi)
1056 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1061 if (molecule_str && *molecule_str &&
1062 strcmp(molecule_str, "(default)")) /* try external PDB files */
1064 /* The -molecule option can point to a .pdb file, or to
1065 a directory of them.
1067 int wire = MI_IS_WIREFRAME(mi);
1074 if (!stat (molecule_str, &st) &&
1075 S_ISDIR (st.st_mode))
1079 struct dirent *dentry;
1081 pdb_dir = opendir (molecule_str);
1084 sprintf (buf, "%.100s: %.100s", progname, molecule_str);
1090 fprintf (stderr, "%s: directory %s\n", progname, molecule_str);
1094 files = (char **) calloc (sizeof(*files), list_size);
1096 while ((dentry = readdir (pdb_dir)))
1098 int L = strlen (dentry->d_name);
1099 if (L > 4 && !strcasecmp (dentry->d_name + L - 4, ".pdb"))
1102 if (nfiles >= list_size-1)
1104 list_size = (list_size + 10) * 1.2;
1106 realloc (files, list_size * sizeof(*files));
1110 fprintf (stderr, "%s: out of memory (%d files)\n",
1116 fn = (char *) malloc (strlen (molecule_str) + L + 10);
1118 strcpy (fn, molecule_str);
1119 if (fn[strlen(fn)-1] != '/') strcat (fn, "/");
1120 strcat (fn, dentry->d_name);
1121 files[nfiles++] = fn;
1123 fprintf (stderr, "%s: file %s\n", progname, fn);
1129 fprintf (stderr, "%s: no .pdb files in directory %s\n",
1130 progname, molecule_str);
1134 files = (char **) malloc (sizeof (*files));
1136 files[0] = strdup (molecule_str);
1138 fprintf (stderr, "%s: file %s\n", progname, molecule_str);
1141 mc->nmolecules = nfiles;
1142 mc->molecules = (molecule *) calloc (sizeof (molecule), mc->nmolecules);
1144 for (i = 0; i < mc->nmolecules; i++)
1147 fprintf (stderr, "%s: reading %s\n", progname, files[i]);
1148 if (!parse_pdb_file (&mc->molecules[molecule_ctr], files[i]))
1150 if ((wire || !do_atoms) &&
1152 mc->molecules[molecule_ctr].nbonds == 0)
1154 /* If we're not drawing atoms (e.g., wireframe mode), and
1155 there is no bond info, then make sure labels are turned
1156 on, or we'll be looking at a black screen... */
1157 fprintf (stderr, "%s: %s: no bonds: turning -label on.\n",
1158 progname, files[i]);
1169 mc->nmolecules = molecule_ctr;
1171 # endif /* LOAD_FILES */
1173 if (mc->nmolecules == 0) /* do the builtins if no files */
1175 mc->nmolecules = countof(builtin_pdb_data);
1176 mc->molecules = (molecule *) calloc (sizeof (molecule), mc->nmolecules);
1177 for (i = 0; i < mc->nmolecules; i++)
1180 sprintf (name, "<builtin-%d>", i);
1181 if (verbose_p) fprintf (stderr, "%s: reading %s\n", progname, name);
1182 parse_pdb_data (&mc->molecules[i], builtin_pdb_data[i], name, 1);
1186 for (i = 0; i < mc->nmolecules; i++)
1188 generate_molecule_formula (&mc->molecules[i]);
1189 insert_vertical_whitespace ((char *) mc->molecules[i].label);
1195 /* Window management, etc
1198 reshape_molecule (ModeInfo *mi, int width, int height)
1200 GLfloat h = (GLfloat) height / (GLfloat) width;
1202 glViewport (0, 0, (GLint) width, (GLint) height);
1204 glMatrixMode(GL_PROJECTION);
1206 gluPerspective (30.0, 1/h, 20.0, 100.0);
1208 glMatrixMode(GL_MODELVIEW);
1210 gluLookAt( 0.0, 0.0, 30.0,
1214 glClear(GL_COLOR_BUFFER_BIT);
1219 gl_init (ModeInfo *mi)
1221 static const GLfloat pos[4] = {1.0, 0.4, 0.9, 0.0};
1222 static const GLfloat amb[4] = {0.0, 0.0, 0.0, 1.0};
1223 static const GLfloat dif[4] = {0.8, 0.8, 0.8, 1.0};
1224 static const GLfloat spc[4] = {1.0, 1.0, 1.0, 1.0};
1225 glLightfv(GL_LIGHT0, GL_POSITION, pos);
1226 glLightfv(GL_LIGHT0, GL_AMBIENT, amb);
1227 glLightfv(GL_LIGHT0, GL_DIFFUSE, dif);
1228 glLightfv(GL_LIGHT0, GL_SPECULAR, spc);
1233 startup_blurb (ModeInfo *mi)
1235 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1236 const char *s = "Constructing molecules...";
1237 print_gl_string (mi->dpy,
1238 # ifdef HAVE_GLBITMAP
1239 mc->xfont3, mc->font3_dlist,
1243 mi->xgwa.width, mi->xgwa.height,
1244 10, mi->xgwa.height - 10,
1247 glXSwapBuffers(MI_DISPLAY(mi), MI_WINDOW(mi));
1251 molecule_handle_event (ModeInfo *mi, XEvent *event)
1253 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1255 if (event->xany.type == ButtonPress &&
1256 event->xbutton.button == Button1)
1258 mc->button_down_p = True;
1259 gltrackball_start (mc->trackball,
1260 event->xbutton.x, event->xbutton.y,
1261 MI_WIDTH (mi), MI_HEIGHT (mi));
1264 else if (event->xany.type == ButtonRelease &&
1265 event->xbutton.button == Button1)
1267 mc->button_down_p = False;
1270 else if (event->xany.type == ButtonPress &&
1271 (event->xbutton.button == Button4 ||
1272 event->xbutton.button == Button5 ||
1273 event->xbutton.button == Button6 ||
1274 event->xbutton.button == Button7))
1276 gltrackball_mousewheel (mc->trackball, event->xbutton.button, 10,
1277 !!event->xbutton.state);
1280 else if (event->xany.type == KeyPress)
1284 XLookupString (&event->xkey, &c, 1, &keysym, 0);
1286 if (c == ' ' || c == '\t' || c == '\r' || c == '\n')
1288 GLfloat speed = 4.0;
1290 mc->mode_tick = 10 * speed;
1294 else if (event->xany.type == MotionNotify &&
1297 gltrackball_track (mc->trackball,
1298 event->xmotion.x, event->xmotion.y,
1299 MI_WIDTH (mi), MI_HEIGHT (mi));
1308 init_molecule (ModeInfo *mi)
1310 molecule_configuration *mc;
1314 mcs = (molecule_configuration *)
1315 calloc (MI_NUM_SCREENS(mi), sizeof (molecule_configuration));
1317 fprintf(stderr, "%s: out of memory\n", progname);
1322 mc = &mcs[MI_SCREEN(mi)];
1324 if ((mc->glx_context = init_GL(mi)) != NULL) {
1326 reshape_molecule (mi, MI_WIDTH(mi), MI_HEIGHT(mi));
1332 wire = MI_IS_WIREFRAME(mi);
1335 Bool spinx=False, spiny=False, spinz=False;
1336 double spin_speed = 0.5;
1337 double spin_accel = 0.3;
1338 double wander_speed = 0.01;
1343 if (*s == 'x' || *s == 'X') spinx = True;
1344 else if (*s == 'y' || *s == 'Y') spiny = True;
1345 else if (*s == 'z' || *s == 'Z') spinz = True;
1346 else if (*s == '0') ;
1350 "%s: spin must contain only the characters X, Y, or Z (not \"%s\")\n",
1357 mc->rot = make_rotator (spinx ? spin_speed : 0,
1358 spiny ? spin_speed : 0,
1359 spinz ? spin_speed : 0,
1361 do_wander ? wander_speed : 0,
1362 (spinx && spiny && spinz));
1363 mc->trackball = gltrackball_init ();
1366 orig_do_labels = do_labels;
1367 orig_do_atoms = do_atoms;
1368 orig_do_bonds = do_bonds;
1369 orig_do_shells = do_shells;
1370 orig_wire = MI_IS_WIREFRAME(mi);
1372 mc->molecule_dlist = glGenLists(1);
1374 mc->shell_dlist = glGenLists(1);
1376 load_molecules (mi);
1377 mc->which = random() % mc->nmolecules;
1379 mc->no_label_threshold = get_float_resource (mi->dpy, "noLabelThreshold",
1380 "NoLabelThreshold");
1381 mc->wireframe_threshold = get_float_resource (mi->dpy, "wireframeThreshold",
1382 "WireframeThreshold");
1390 /* Put the labels on the atoms.
1391 This can't be a part of the display list because of the games
1392 we play with the translation matrix.
1395 draw_labels (ModeInfo *mi)
1397 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1398 int wire = MI_IS_WIREFRAME(mi);
1399 molecule *m = &mc->molecules[mc->which];
1400 # ifdef HAVE_GLBITMAP
1401 XFontStruct *xfont = (mc->scale_down ? mc->xfont2 : mc->xfont1);
1402 GLuint font_dlist = (mc->scale_down ? mc->font2_dlist : mc->font1_dlist);
1404 texture_font_data *font_data = mc->font1_data; /* don't scale down */
1412 glDisable (GL_LIGHTING); /* don't light fonts */
1414 for (i = 0; i < m->natoms; i++)
1416 molecule_atom *a = &m->atoms[i];
1417 GLfloat size = atom_size (a);
1423 set_atom_color (mi, a, True, 1);
1425 /* First, we translate the origin to the center of the atom.
1427 Then we retrieve the prevailing modelview matrix (which
1428 includes any rotation, wandering, and user-trackball-rolling
1431 We set the top 3x3 cells of that matrix to be the identity
1432 matrix. This removes all rotation from the matrix, while
1433 leaving the translation alone. This has the effect of
1434 leaving the prevailing coordinate system perpendicular to
1435 the camera view: were we to draw a square face, it would
1436 be in the plane of the screen.
1438 Now we translate by `size' toward the viewer -- so that the
1439 origin is *just in front* of the ball.
1441 Then we draw the label text, allowing the depth buffer to
1442 do its work: that way, labels on atoms will be occluded
1443 properly when other atoms move in front of them.
1445 This technique (of neutralizing rotation relative to the
1446 observer, after both rotations and translations have been
1447 applied) is known as "billboarding".
1450 glTranslatef(a->x, a->y, a->z); /* get matrix */
1451 glGetFloatv (GL_MODELVIEW_MATRIX, &m[0][0]); /* load rot. identity */
1452 m[0][0] = 1; m[1][0] = 0; m[2][0] = 0;
1453 m[0][1] = 0; m[1][1] = 1; m[2][1] = 0;
1454 m[0][2] = 0; m[1][2] = 0; m[2][2] = 1;
1455 glLoadIdentity(); /* reset modelview */
1456 glMultMatrixf (&m[0][0]); /* replace with ours */
1458 glTranslatef (0, 0, (size * 1.1)); /* move toward camera */
1460 glRotatef (current_device_rotation(), 0, 0, 1); /* right side up */
1462 # ifdef HAVE_GLBITMAP
1463 glRasterPos3f (0, 0, 0); /* draw text here */
1465 /* Before drawing the string, shift the origin to center
1466 the text over the origin of the sphere. */
1467 glBitmap (0, 0, 0, 0,
1468 -string_width (xfont, a->label, 0) / 2,
1473 for (j = 0; j < strlen(a->label); j++)
1474 glCallList (font_dlist + (int)(a->label[j]));
1479 int w = texture_string_width (font_data, a->label, &h);
1480 GLfloat s = 1.0 / h;
1481 GLfloat max = 18; /* max point size to avoid pixellated text */
1482 if (h > max) s *= max/h;
1484 glTranslatef (-w/2, h*2/3, 0);
1485 print_gl_string (mi->dpy,
1486 # ifdef HAVE_GLBITMAP
1499 /* More efficient to always call glEnable() with correct values
1500 than to call glPushAttrib()/glPopAttrib(), since reading
1501 attributes from GL does a round-trip and stalls the pipeline.
1504 glEnable (GL_LIGHTING);
1509 pick_new_molecule (ModeInfo *mi, time_t last)
1511 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1513 if (mc->nmolecules == 1)
1515 if (last != 0) return;
1520 mc->which = random() % mc->nmolecules;
1525 while (n == mc->which)
1526 n = random() % mc->nmolecules;
1532 char *name = strdup (mc->molecules[mc->which].label);
1533 char *s = strpbrk (name, "\r\n");
1535 fprintf (stderr, "%s: drawing %s (%d)\n", progname, name, mc->which);
1539 mc->polygon_count = 0;
1541 glNewList (mc->molecule_dlist, GL_COMPILE);
1542 ensure_bounding_box_visible (mi);
1544 do_labels = orig_do_labels;
1545 do_atoms = orig_do_atoms;
1546 do_bonds = orig_do_bonds;
1547 do_shells = orig_do_shells;
1548 MI_IS_WIREFRAME(mi) = orig_wire;
1550 if (mc->molecule_size > mc->no_label_threshold)
1552 if (mc->molecule_size > mc->wireframe_threshold)
1553 MI_IS_WIREFRAME(mi) = 1;
1555 if (MI_IS_WIREFRAME(mi))
1556 do_bonds = 1, do_shells = 0;
1561 if (! (do_bonds || do_atoms || do_labels))
1563 /* Make sure *something* shows up! */
1564 MI_IS_WIREFRAME(mi) = 1;
1568 build_molecule (mi, False);
1573 glNewList (mc->shell_dlist, GL_COMPILE);
1574 ensure_bounding_box_visible (mi);
1580 build_molecule (mi, True);
1583 do_bonds = orig_do_bonds;
1584 do_atoms = orig_do_atoms;
1585 do_labels = orig_do_labels;
1591 draw_molecule (ModeInfo *mi)
1593 time_t now = time ((time_t *) 0);
1594 GLfloat speed = 4.0; /* speed at which the zoom out/in happens */
1596 molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1597 Display *dpy = MI_DISPLAY(mi);
1598 Window window = MI_WINDOW(mi);
1600 if (!mc->glx_context)
1603 glXMakeCurrent(MI_DISPLAY(mi), MI_WINDOW(mi), *(mc->glx_context));
1605 if (mc->draw_time == 0)
1607 pick_new_molecule (mi, mc->draw_time);
1608 mc->draw_time = now;
1610 else if (mc->mode == 0)
1612 if (mc->draw_tick++ > 10)
1614 time_t now = time((time_t *) 0);
1615 if (mc->draw_time == 0) mc->draw_time = now;
1618 if (!mc->button_down_p &&
1619 mc->nmolecules > 1 &&
1620 mc->draw_time + timeout <= now)
1622 /* randomize molecules every -timeout seconds */
1623 mc->mode = 1; /* go out */
1624 mc->mode_tick = 10 * speed;
1625 mc->draw_time = now;
1629 else if (mc->mode == 1) /* out */
1631 if (--mc->mode_tick <= 0)
1633 mc->mode_tick = 10 * speed;
1634 mc->mode = 2; /* go in */
1635 pick_new_molecule (mi, mc->draw_time);
1636 mc->draw_time = now;
1639 else if (mc->mode == 2) /* in */
1641 if (--mc->mode_tick <= 0)
1642 mc->mode = 0; /* normal */
1648 glScalef(1.1, 1.1, 1.1);
1652 get_position (mc->rot, &x, &y, &z, !mc->button_down_p);
1653 glTranslatef((x - 0.5) * 9,
1657 /* Do it twice because we don't track the device's orientation. */
1658 glRotatef( current_device_rotation(), 0, 0, 1);
1659 gltrackball_rotate (mc->trackball);
1660 glRotatef(-current_device_rotation(), 0, 0, 1);
1662 get_rotation (mc->rot, &x, &y, &z, !mc->button_down_p);
1663 glRotatef (x * 360, 1.0, 0.0, 0.0);
1664 glRotatef (y * 360, 0.0, 1.0, 0.0);
1665 glRotatef (z * 360, 0.0, 0.0, 1.0);
1668 glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
1672 GLfloat s = (mc->mode == 1
1673 ? mc->mode_tick / (10 * speed)
1674 : ((10 * speed) - mc->mode_tick + 1) / (10 * speed));
1679 glCallList (mc->molecule_dlist);
1683 molecule *m = &mc->molecules[mc->which];
1687 /* This can't go in the display list, or the characters are spaced
1688 wrongly when the window is resized. */
1689 if (do_titles && m->label && *m->label)
1691 set_atom_color (mi, 0, True, 1);
1692 print_gl_string (mi->dpy,
1693 # ifdef HAVE_GLBITMAP
1694 mc->xfont3, mc->font3_dlist,
1698 mi->xgwa.width, mi->xgwa.height,
1699 10, mi->xgwa.height - 10,
1707 glColorMask (GL_FALSE, GL_FALSE, GL_FALSE, GL_FALSE);
1709 glCallList (mc->shell_dlist);
1711 glColorMask (GL_TRUE, GL_TRUE, GL_TRUE, GL_TRUE);
1713 glDepthFunc (GL_EQUAL);
1714 glEnable (GL_BLEND);
1715 glBlendFunc (GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
1717 glCallList (mc->shell_dlist);
1719 glDepthFunc (GL_LESS);
1720 glDisable (GL_BLEND);
1725 mi->polygon_count = mc->polygon_count;
1727 if (mi->fps_p) do_fps (mi);
1730 glXSwapBuffers(dpy, window);
1733 XSCREENSAVER_MODULE ("Molecule", molecule)