ftp://ftp.krokus.ru/pub/OpenBSD/distfiles/xscreensaver-4.22.tar.gz
[xscreensaver] / hacks / glx / molecule.c
1 /* molecule, Copyright (c) 2001-2005 Jamie Zawinski <jwz@jwz.org>
2  * Draws molecules, based on coordinates from PDB (Protein Data Base) files.
3  *
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 
10  * implied warranty.
11  */
12
13
14 /* Documentation on the PDB file format:
15    http://www.rcsb.org/pdb/docs/format/pdbguide2.2/guide2.2_frame.html
16
17    Good source of PDB files:
18    http://www.sci.ouc.bc.ca/chem/molecule/molecule.html
19  */
20
21 #include <sys/types.h>
22 #include <sys/stat.h>
23 #include <unistd.h>
24 #include <dirent.h>
25 #include <X11/Intrinsic.h>
26
27 #define PROGCLASS       "Molecule"
28 #define HACK_INIT       init_molecule
29 #define HACK_DRAW       draw_molecule
30 #define HACK_RESHAPE    reshape_molecule
31 #define HACK_HANDLE_EVENT molecule_handle_event
32 #define EVENT_MASK      PointerMotionMask
33 #define molecule_opts   xlockmore_opts
34
35 #define DEF_TIMEOUT     "20"
36 #define DEF_SPIN        "XYZ"
37 #define DEF_WANDER      "False"
38 #define DEF_LABELS      "True"
39 #define DEF_TITLES      "True"
40 #define DEF_ATOMS       "True"
41 #define DEF_BONDS       "True"
42 #define DEF_SHELLS      "False"
43 #define DEF_BBOX        "False"
44 #define DEF_SHELL_ALPHA "0.3"
45 #define DEF_MOLECULE    "(default)"
46 #define DEF_VERBOSE     "False"
47
48 #define DEFAULTS        "*delay:        10000         \n" \
49                         "*showFPS:      False         \n" \
50                         "*wireframe:    False         \n" \
51                         "*atomFont:   -*-times-bold-r-normal-*-240-*\n" \
52                         "*titleFont:  -*-times-bold-r-normal-*-180-*\n" \
53                         "*noLabelThreshold:    30     \n" \
54                         "*wireframeThreshold:  150    \n" \
55
56
57 #undef countof
58 #define countof(x) (sizeof((x))/sizeof((*x)))
59
60 #include "xlockmore.h"
61 #include "colors.h"
62 #include "sphere.h"
63 #include "tube.h"
64 #include "glxfonts.h"
65 #include "rotator.h"
66 #include "gltrackball.h"
67
68 #ifdef USE_GL /* whole file */
69
70 #include <stdlib.h>
71 #include <ctype.h>
72 #include <time.h>
73 #include <sys/time.h>
74 #include <GL/glu.h>
75
76 #define SPHERE_SLICES 24  /* how densely to render spheres */
77 #define SPHERE_STACKS 12
78
79 #define SMOOTH_TUBE       /* whether to have smooth or faceted tubes */
80
81 #ifdef SMOOTH_TUBE
82 # define TUBE_FACES  12   /* how densely to render tubes */
83 #else
84 # define TUBE_FACES  8
85 #endif
86
87 static int scale_down;
88 #define SPHERE_SLICES_2  7
89 #define SPHERE_STACKS_2  4
90 #define TUBE_FACES_2     3
91
92
93 # ifdef __GNUC__
94   __extension__  /* don't warn about "string length is greater than the length
95                     ISO C89 compilers are required to support" when includng
96                     the following data file... */
97 # endif
98 const char * const builtin_pdb_data[] = {
99 # include "molecules.h"
100 };
101
102
103 typedef struct {
104   const char *name;
105   GLfloat size, size2;
106   const char *color;
107   const char *text_color;
108   GLfloat gl_color[8];
109 } atom_data;
110
111
112 /* These are the traditional colors used to render these atoms,
113    and their approximate size in angstroms.
114  */
115 static atom_data all_atom_data[] = {
116   { "H",    1.17,  0,  "White",           "Grey70",        { 0, }},
117   { "C",    1.75,  0,  "Grey60",          "White",         { 0, }},
118   { "CA",   1.80,  0,  "Blue",            "LightBlue",     { 0, }},
119   { "N",    1.55,  0,  "LightSteelBlue3", "SlateBlue1",    { 0, }},
120   { "O",    1.40,  0,  "Red",             "LightPink",     { 0, }},
121   { "P",    1.28,  0,  "MediumPurple",    "PaleVioletRed", { 0, }},
122   { "S",    1.80,  0,  "Yellow4",         "Yellow1",       { 0, }},
123   { "bond", 0,     0,  "Grey70",          "Yellow1",       { 0, }},
124   { "*",    1.40,  0,  "Green4",          "LightGreen",    { 0, }}
125 };
126
127
128 typedef struct {
129   int id;               /* sequence number in the PDB file */
130   const char *label;    /* The atom name */
131   GLfloat x, y, z;      /* position in 3-space (angstroms) */
132   atom_data *data;      /* computed: which style of atom this is */
133 } molecule_atom;
134
135 typedef struct {
136   int from, to;         /* atom sequence numbers */
137   int strength;         /* how many bonds are between these two atoms */
138 } molecule_bond;
139
140
141 typedef struct {
142   const char *label;            /* description of this compound */
143   int natoms, atoms_size;
144   int nbonds, bonds_size;
145   molecule_atom *atoms;
146   molecule_bond *bonds;
147 } molecule;
148
149
150 typedef struct {
151   GLXContext *glx_context;
152   rotator *rot;
153   trackball_state *trackball;
154   Bool button_down_p;
155
156   GLfloat molecule_size;           /* max dimension of molecule bounding box */
157
158   GLfloat no_label_threshold;      /* Things happen when molecules are huge */
159   GLfloat wireframe_threshold;
160
161   int which;                       /* which of the molecules is being shown */
162   int nmolecules;
163   molecule *molecules;
164
165   int mode;  /* 0 = normal, 1 = out, 2 = in */
166   int mode_tick;
167
168   GLuint molecule_dlist;
169   GLuint shell_dlist;
170
171   XFontStruct *xfont1, *xfont2;
172   GLuint font1_dlist, font2_dlist;
173   int polygon_count;
174
175 } molecule_configuration;
176
177
178 static molecule_configuration *mcs = NULL;
179
180 static int timeout;
181 static char *molecule_str;
182 static char *do_spin;
183 static Bool do_wander;
184 static Bool do_titles;
185 static Bool do_labels;
186 static Bool do_atoms;
187 static Bool do_bonds;
188 static Bool do_shells;
189 static Bool do_bbox;
190 static Bool verbose_p;
191 static GLfloat shell_alpha;
192
193 /* saved to reset */
194 static Bool orig_do_labels, orig_do_atoms, orig_do_bonds, orig_do_shells,
195     orig_wire;
196
197
198 static XrmOptionDescRec opts[] = {
199   { "-molecule",        ".molecule",    XrmoptionSepArg, 0 },
200   { "-timeout",         ".timeout",     XrmoptionSepArg, 0 },
201   { "-spin",            ".spin",        XrmoptionSepArg, 0 },
202   { "+spin",            ".spin",        XrmoptionNoArg, "" },
203   { "-wander",          ".wander",      XrmoptionNoArg, "True" },
204   { "+wander",          ".wander",      XrmoptionNoArg, "False" },
205   { "-labels",          ".labels",      XrmoptionNoArg, "True" },
206   { "+labels",          ".labels",      XrmoptionNoArg, "False" },
207   { "-titles",          ".titles",      XrmoptionNoArg, "True" },
208   { "+titles",          ".titles",      XrmoptionNoArg, "False" },
209   { "-atoms",           ".atoms",       XrmoptionNoArg, "True" },
210   { "+atoms",           ".atoms",       XrmoptionNoArg, "False" },
211   { "-bonds",           ".bonds",       XrmoptionNoArg, "True" },
212   { "+bonds",           ".bonds",       XrmoptionNoArg, "False" },
213   { "-shells",          ".eshells",     XrmoptionNoArg, "True" },
214   { "+shells",          ".eshells",     XrmoptionNoArg, "False" },
215   { "-shell-alpha",     ".shellAlpha",  XrmoptionSepArg, 0 },
216   { "-bbox",            ".bbox",        XrmoptionNoArg, "True" },
217   { "+bbox",            ".bbox",        XrmoptionNoArg, "False" },
218   { "-verbose",         ".verbose",     XrmoptionNoArg, "True" },
219 };
220
221 static argtype vars[] = {
222   {&molecule_str, "molecule",   "Molecule",     DEF_MOLECULE, t_String},
223   {&timeout,      "timeout",    "Seconds",      DEF_TIMEOUT,  t_Int},
224   {&do_spin,      "spin",       "Spin",         DEF_SPIN,     t_String},
225   {&do_wander,    "wander",     "Wander",       DEF_WANDER,   t_Bool},
226   {&do_atoms,     "atoms",      "Atoms",        DEF_ATOMS,    t_Bool},
227   {&do_bonds,     "bonds",      "Bonds",        DEF_BONDS,    t_Bool},
228   {&do_shells,    "eshells",    "EShells",      DEF_SHELLS,   t_Bool},
229   {&do_labels,    "labels",     "Labels",       DEF_LABELS,   t_Bool},
230   {&do_titles,    "titles",     "Titles",       DEF_TITLES,   t_Bool},
231   {&do_bbox,      "bbox",       "BBox",         DEF_BBOX,     t_Bool},
232   {&shell_alpha,  "shellAlpha", "ShellAlpha",   DEF_SHELL_ALPHA, t_Float},
233   {&verbose_p,    "verbose",    "Verbose",      DEF_VERBOSE,  t_Bool},
234 };
235
236 ModeSpecOpt molecule_opts = {countof(opts), opts, countof(vars), vars, NULL};
237
238
239
240 \f
241 /* shapes */
242
243 static int
244 sphere (GLfloat x, GLfloat y, GLfloat z, GLfloat diameter, Bool wire)
245 {
246   int stacks = (scale_down ? SPHERE_STACKS_2 : SPHERE_STACKS);
247   int slices = (scale_down ? SPHERE_SLICES_2 : SPHERE_SLICES);
248
249   glPushMatrix ();
250   glTranslatef (x, y, z);
251   glScalef (diameter, diameter, diameter);
252   unit_sphere (stacks, slices, wire);
253   glPopMatrix ();
254
255   return stacks * slices;
256 }
257
258
259 static void
260 load_fonts (ModeInfo *mi)
261 {
262   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
263   load_font (mi->dpy, "atomFont",  &mc->xfont1, &mc->font1_dlist);
264   load_font (mi->dpy, "titleFont", &mc->xfont2, &mc->font2_dlist);
265 }
266
267
268 static atom_data *
269 get_atom_data (const char *atom_name)
270 {
271   int i;
272   atom_data *d = 0;
273   char *n = strdup (atom_name);
274   char *n2 = n;
275   int L;
276
277   while (!isalpha(*n)) n++;
278   L = strlen(n);
279   while (L > 0 && !isalpha(n[L-1]))
280     n[--L] = 0;
281
282   for (i = 0; i < countof(all_atom_data); i++)
283     {
284       d = &all_atom_data[i];
285       if (!strcasecmp (n, all_atom_data[i].name))
286         break;
287     }
288
289   free (n2);
290   return d;
291 }
292
293
294 static void
295 set_atom_color (ModeInfo *mi, molecule_atom *a, Bool font_p, GLfloat alpha)
296 {
297   atom_data *d;
298   GLfloat *gl_color;
299
300   if (a)
301     d = a->data;
302   else
303     {
304       static atom_data *def_data = 0;
305       if (!def_data) def_data = get_atom_data ("bond");
306       d = def_data;
307     }
308
309   gl_color = (!font_p ? d->gl_color : (d->gl_color + 4));
310
311   if (gl_color[3] == 0)
312     {
313       const char *string = !font_p ? d->color : d->text_color;
314       XColor xcolor;
315       if (!XParseColor (mi->dpy, mi->xgwa.colormap, string, &xcolor))
316         {
317           fprintf (stderr, "%s: unparsable color in %s: %s\n", progname,
318                    (a ? a->label : d->name), string);
319           exit (1);
320         }
321
322       gl_color[0] = xcolor.red   / 65536.0;
323       gl_color[1] = xcolor.green / 65536.0;
324       gl_color[2] = xcolor.blue  / 65536.0;
325     }
326   
327   gl_color[3] = alpha;
328
329   if (font_p)
330     glColor4f (gl_color[0], gl_color[1], gl_color[2], gl_color[3]);
331   else
332     glMaterialfv (GL_FRONT, GL_AMBIENT_AND_DIFFUSE, gl_color);
333 }
334
335
336 static GLfloat
337 atom_size (molecule_atom *a)
338 {
339   if (do_bonds)
340     {
341       if (a->data->size2 == 0)
342         {
343           /* let the molecules have the same relative sizes, but scale
344              them to a smaller range, so that the bond-tubes are
345              actually visible...
346            */
347           GLfloat bot = 0.4;
348           GLfloat top = 0.6;
349           GLfloat min = 1.17;
350           GLfloat max = 1.80;
351           GLfloat ratio = (a->data->size - min) / (max - min);
352           a->data->size2 = bot + (ratio * (top - bot));
353         }
354       return a->data->size2;
355     }
356   else
357     return a->data->size;
358 }
359
360
361 static molecule_atom *
362 get_atom (molecule_atom *atoms, int natoms, int id)
363 {
364   int i;
365
366   /* quick short-circuit */
367   if (id < natoms)
368     {
369       if (atoms[id].id == id)
370         return &atoms[id];
371       if (id > 0 && atoms[id-1].id == id)
372         return &atoms[id-1];
373       if (id < natoms-1 && atoms[id+1].id == id)
374         return &atoms[id+1];
375     }
376
377   for (i = 0; i < natoms; i++)
378     if (id == atoms[i].id)
379       return &atoms[i];
380
381   fprintf (stderr, "%s: no atom %d\n", progname, id);
382   abort();
383 }
384
385
386 static void
387 molecule_bounding_box (ModeInfo *mi,
388                        GLfloat *x1, GLfloat *y1, GLfloat *z1,
389                        GLfloat *x2, GLfloat *y2, GLfloat *z2)
390 {
391   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
392   molecule *m = &mc->molecules[mc->which];
393   int i;
394
395   if (m->natoms == 0)
396     {
397       *x1 = *y1 = *z1 = *x2 = *y2 = *z2 = 0;
398     }
399   else
400     {
401       *x1 = *x2 = m->atoms[0].x;
402       *y1 = *y2 = m->atoms[0].y;
403       *z1 = *z2 = m->atoms[0].z;
404     }
405
406   for (i = 1; i < m->natoms; i++)
407     {
408       if (m->atoms[i].x < *x1) *x1 = m->atoms[i].x;
409       if (m->atoms[i].y < *y1) *y1 = m->atoms[i].y;
410       if (m->atoms[i].z < *z1) *z1 = m->atoms[i].z;
411
412       if (m->atoms[i].x > *x2) *x2 = m->atoms[i].x;
413       if (m->atoms[i].y > *y2) *y2 = m->atoms[i].y;
414       if (m->atoms[i].z > *z2) *z2 = m->atoms[i].z;
415     }
416
417   *x1 -= 1.5;
418   *y1 -= 1.5;
419   *z1 -= 1.5;
420   *x2 += 1.5;
421   *y2 += 1.5;
422   *z2 += 1.5;
423 }
424
425
426 static void
427 draw_bounding_box (ModeInfo *mi)
428 {
429   static GLfloat c1[4] = { 0.2, 0.2, 0.4, 1.0 };
430   static GLfloat c2[4] = { 1.0, 0.0, 0.0, 1.0 };
431   int wire = MI_IS_WIREFRAME(mi);
432   GLfloat x1, y1, z1, x2, y2, z2;
433   molecule_bounding_box (mi, &x1, &y1, &z1, &x2, &y2, &z2);
434
435   glMaterialfv (GL_FRONT, GL_AMBIENT_AND_DIFFUSE, c1);
436   glFrontFace(GL_CCW);
437
438   glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
439   glNormal3f(0, 1, 0);
440   glVertex3f(x1, y1, z1); glVertex3f(x1, y1, z2);
441   glVertex3f(x2, y1, z2); glVertex3f(x2, y1, z1);
442   glEnd();
443   glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
444   glNormal3f(0, -1, 0);
445   glVertex3f(x2, y2, z1); glVertex3f(x2, y2, z2);
446   glVertex3f(x1, y2, z2); glVertex3f(x1, y2, z1);
447   glEnd();
448   glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
449   glNormal3f(0, 0, 1);
450   glVertex3f(x1, y1, z1); glVertex3f(x2, y1, z1);
451   glVertex3f(x2, y2, z1); glVertex3f(x1, y2, z1);
452   glEnd();
453   glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
454   glNormal3f(0, 0, -1);
455   glVertex3f(x1, y2, z2); glVertex3f(x2, y2, z2);
456   glVertex3f(x2, y1, z2); glVertex3f(x1, y1, z2);
457   glEnd();
458   glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
459   glNormal3f(1, 0, 0);
460   glVertex3f(x1, y2, z1); glVertex3f(x1, y2, z2);
461   glVertex3f(x1, y1, z2); glVertex3f(x1, y1, z1);
462   glEnd();
463   glBegin(wire ? GL_LINE_LOOP : GL_QUADS);
464   glNormal3f(-1, 0, 0);
465   glVertex3f(x2, y1, z1); glVertex3f(x2, y1, z2);
466   glVertex3f(x2, y2, z2); glVertex3f(x2, y2, z1);
467   glEnd();
468
469   glPushAttrib (GL_LIGHTING);
470   glDisable (GL_LIGHTING);
471
472   glColor3f (c2[0], c2[1], c2[2]);
473   glBegin(GL_LINES);
474   if (x1 > 0) x1 = 0; if (x2 < 0) x2 = 0;
475   if (y1 > 0) y1 = 0; if (y2 < 0) y2 = 0;
476   if (z1 > 0) z1 = 0; if (z2 < 0) z2 = 0;
477   glVertex3f(x1, 0,  0);  glVertex3f(x2, 0,  0); 
478   glVertex3f(0 , y1, 0);  glVertex3f(0,  y2, 0); 
479   glVertex3f(0,  0,  z1); glVertex3f(0,  0,  z2); 
480   glEnd();
481
482   glPopAttrib();
483 }
484
485
486 /* Since PDB files don't always have the molecule centered around the
487    origin, and since some molecules are pretty large, scale and/or
488    translate so that the whole molecule is visible in the window.
489  */
490 static void
491 ensure_bounding_box_visible (ModeInfo *mi)
492 {
493   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
494
495   GLfloat x1, y1, z1, x2, y2, z2;
496   GLfloat w, h, d;
497   GLfloat size;
498   GLfloat max_size = 10;  /* don't bother scaling down if the molecule
499                              is already smaller than this */
500
501   molecule_bounding_box (mi, &x1, &y1, &z1, &x2, &y2, &z2);
502   w = x2-x1;
503   h = y2-y1;
504   d = z2-z1;
505
506   size = (w > h ? w : h);
507   size = (size > d ? size : d);
508
509   mc->molecule_size = size;
510
511   scale_down = 0;
512
513   if (size > max_size)
514     {
515       GLfloat scale = max_size / size;
516       glScalef (scale, scale, scale);
517
518       scale_down = scale < 0.3;
519     }
520
521   glTranslatef (-(x1 + w/2),
522                 -(y1 + h/2),
523                 -(z1 + d/2));
524 }
525
526
527
528 /* Constructs the GL shapes of the current molecule
529  */
530 static void
531 build_molecule (ModeInfo *mi, Bool transparent_p)
532 {
533   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
534   int wire = MI_IS_WIREFRAME(mi);
535   int i;
536   GLfloat alpha = transparent_p ? shell_alpha : 1.0;
537   int polys = 0;
538
539   molecule *m = &mc->molecules[mc->which];
540
541   if (wire)
542     {
543       glDisable(GL_CULL_FACE);
544       glDisable(GL_LIGHTING);
545       glDisable(GL_LIGHT0);
546       glDisable(GL_DEPTH_TEST);
547       glDisable(GL_NORMALIZE);
548       glDisable(GL_CULL_FACE);
549     }
550   else
551     {
552       glEnable(GL_CULL_FACE);
553       glEnable(GL_LIGHTING);
554       glEnable(GL_LIGHT0);
555       glEnable(GL_DEPTH_TEST);
556       glEnable(GL_NORMALIZE);
557       glEnable(GL_CULL_FACE);
558     }
559
560   if (!wire)
561     set_atom_color (mi, 0, False, alpha);
562
563   if (do_bonds)
564     for (i = 0; i < m->nbonds; i++)
565       {
566         molecule_bond *b = &m->bonds[i];
567         molecule_atom *from = get_atom (m->atoms, m->natoms, b->from);
568         molecule_atom *to   = get_atom (m->atoms, m->natoms, b->to);
569
570         if (wire)
571           {
572             glBegin(GL_LINES);
573             glVertex3f(from->x, from->y, from->z);
574             glVertex3f(to->x,   to->y,   to->z);
575             glEnd();
576             polys++;
577           }
578         else
579           {
580             int faces = (scale_down ? TUBE_FACES_2 : TUBE_FACES);
581 # ifdef SMOOTH_TUBE
582             int smooth = True;
583 # else
584             int smooth = False;
585 # endif
586             GLfloat thickness = 0.07 * b->strength;
587             GLfloat cap_size = 0.03;
588             if (thickness > 0.3)
589               thickness = 0.3;
590
591             tube (from->x, from->y, from->z,
592                   to->x,   to->y,   to->z,
593                   thickness, cap_size,
594                   faces, smooth, (!do_atoms || do_shells), wire);
595             polys += faces;
596           }
597       }
598
599   if (!wire && do_atoms)
600     for (i = 0; i < m->natoms; i++)
601       {
602         molecule_atom *a = &m->atoms[i];
603         GLfloat size = atom_size (a);
604         set_atom_color (mi, a, False, alpha);
605         polys += sphere (a->x, a->y, a->z, size, wire);
606       }
607
608   if (do_bbox && !transparent_p)
609     {
610       draw_bounding_box (mi);
611       polys += 4;
612     }
613
614   mc->polygon_count += polys;
615 }
616
617
618 \f
619 /* loading */
620
621 static void
622 push_atom (molecule *m,
623            int id, const char *label,
624            GLfloat x, GLfloat y, GLfloat z)
625 {
626   m->natoms++;
627   if (m->atoms_size < m->natoms)
628     {
629       m->atoms_size += 20;
630       m->atoms = (molecule_atom *) realloc (m->atoms,
631                                             m->atoms_size * sizeof(*m->atoms));
632     }
633   m->atoms[m->natoms-1].id = id;
634   m->atoms[m->natoms-1].label = label;
635   m->atoms[m->natoms-1].x = x;
636   m->atoms[m->natoms-1].y = y;
637   m->atoms[m->natoms-1].z = z;
638   m->atoms[m->natoms-1].data = get_atom_data (label);
639 }
640
641
642 static void
643 push_bond (molecule *m, int from, int to)
644 {
645   int i;
646
647   for (i = 0; i < m->nbonds; i++)
648     if ((m->bonds[i].from == from && m->bonds[i].to   == to) ||
649         (m->bonds[i].to   == from && m->bonds[i].from == to))
650       {
651         m->bonds[i].strength++;
652         return;
653       }
654
655   m->nbonds++;
656   if (m->bonds_size < m->nbonds)
657     {
658       m->bonds_size += 20;
659       m->bonds = (molecule_bond *) realloc (m->bonds,
660                                             m->bonds_size * sizeof(*m->bonds));
661     }
662   m->bonds[m->nbonds-1].from = from;
663   m->bonds[m->nbonds-1].to = to;
664   m->bonds[m->nbonds-1].strength = 1;
665 }
666
667
668
669 /* This function is crap.
670  */
671 static void
672 parse_pdb_data (molecule *m, const char *data, const char *filename, int line)
673 {
674   const char *s = data;
675   char *ss;
676   while (*s)
677     {
678       if ((!m->label || !*m->label) &&
679           (!strncmp (s, "HEADER", 6) || !strncmp (s, "COMPND", 6)))
680         {
681           char *name = calloc (1, 100);
682           char *n2 = name;
683           int L = strlen(s);
684           if (L > 99) L = 99;
685
686           strncpy (n2, s, L);
687           n2 += 7;
688           while (isspace(*n2)) n2++;
689
690           ss = strchr (n2, '\n');
691           if (ss) *ss = 0;
692           ss = strchr (n2, '\r');
693           if (ss) *ss = 0;
694
695           ss = n2+strlen(n2)-1;
696           while (isspace(*ss) && ss > n2)
697             *ss-- = 0;
698
699           if (strlen (n2) > 4 &&
700               !strcmp (n2 + strlen(n2) - 4, ".pdb"))
701             n2[strlen(n2)-4] = 0;
702
703           if (m->label) free ((char *) m->label);
704           m->label = strdup (n2);
705           free (name);
706         }
707       else if (!strncmp (s, "TITLE ", 6) ||
708                !strncmp (s, "HEADER", 6) ||
709                !strncmp (s, "COMPND", 6) ||
710                !strncmp (s, "AUTHOR", 6) ||
711                !strncmp (s, "REVDAT", 6) ||
712                !strncmp (s, "SOURCE", 6) ||
713                !strncmp (s, "EXPDTA", 6) ||
714                !strncmp (s, "JRNL  ", 6) ||
715                !strncmp (s, "REMARK", 6) ||
716                !strncmp (s, "SEQRES", 6) ||
717                !strncmp (s, "HET   ", 6) ||
718                !strncmp (s, "FORMUL", 6) ||
719                !strncmp (s, "CRYST1", 6) ||
720                !strncmp (s, "ORIGX1", 6) ||
721                !strncmp (s, "ORIGX2", 6) ||
722                !strncmp (s, "ORIGX3", 6) ||
723                !strncmp (s, "SCALE1", 6) ||
724                !strncmp (s, "SCALE2", 6) ||
725                !strncmp (s, "SCALE3", 6) ||
726                !strncmp (s, "MASTER", 6) ||
727                !strncmp (s, "KEYWDS", 6) ||
728                !strncmp (s, "DBREF ", 6) ||
729                !strncmp (s, "HETNAM", 6) ||
730                !strncmp (s, "HETSYN", 6) ||
731                !strncmp (s, "HELIX ", 6) ||
732                !strncmp (s, "LINK  ", 6) ||
733                !strncmp (s, "MTRIX1", 6) ||
734                !strncmp (s, "MTRIX2", 6) ||
735                !strncmp (s, "MTRIX3", 6) ||
736                !strncmp (s, "SHEET ", 6) ||
737                !strncmp (s, "CISPEP", 6) ||
738                !strncmp (s, "GENERATED BY", 12) ||
739                !strncmp (s, "TER ", 4) ||
740                !strncmp (s, "END ", 4) ||
741                !strncmp (s, "TER\n", 4) ||
742                !strncmp (s, "END\n", 4) ||
743                !strncmp (s, "\n", 1))
744         /* ignored. */
745         ;
746       else if (!strncmp (s, "ATOM   ", 7))
747         {
748           int id;
749           char *name = (char *) calloc (1, 4);
750           GLfloat x = -999, y = -999, z = -999;
751
752           sscanf (s+7, " %d ", &id);
753
754           strncpy (name, s+12, 3);
755           while (isspace(*name)) name++;
756           ss = name + strlen(name)-1;
757           while (isspace(*ss) && ss > name)
758             *ss-- = 0;
759           ss = name + 1;
760           while(*ss)
761           {
762             *ss = tolower(*ss);
763             ss++;
764           }
765           sscanf (s + 32, " %f %f %f ", &x, &y, &z);
766 /*
767           fprintf (stderr, "%s: %s: %d: atom: %d \"%s\" %9.4f %9.4f %9.4f\n",
768                    progname, filename, line,
769                    id, name, x, y, z);
770 */
771           push_atom (m, id, name, x, y, z);
772         }
773       else if (!strncmp (s, "HETATM ", 7))
774         {
775           int id;
776           char *name = (char *) calloc (1, 4);
777           GLfloat x = -999, y = -999, z = -999;
778
779           sscanf (s+7, " %d ", &id);
780
781           strncpy (name, s+12, 3);
782           while (isspace(*name)) name++;
783           ss = name + strlen(name)-1;
784           while (isspace(*ss) && ss > name)
785             *ss-- = 0;
786           sscanf (s + 30, " %f %f %f ", &x, &y, &z);
787 /*
788           fprintf (stderr, "%s: %s: %d: atom: %d \"%s\" %9.4f %9.4f %9.4f\n",
789                    progname, filename, line,
790                    id, name, x, y, z);
791 */
792           push_atom (m, id, name, x, y, z);
793         }
794       else if (!strncmp (s, "CONECT ", 7))
795         {
796           int atoms[11];
797           int i = sscanf (s + 8, " %d %d %d %d %d %d %d %d %d %d %d %d ",
798                           &atoms[0], &atoms[1], &atoms[2], &atoms[3], 
799                           &atoms[4], &atoms[5], &atoms[6], &atoms[7], 
800                           &atoms[8], &atoms[9], &atoms[10], &atoms[11]);
801           int j;
802           for (j = 1; j < i; j++)
803             if (atoms[j] > 0)
804               {
805 /*
806                 fprintf (stderr, "%s: %s: %d: bond: %d %d\n",
807                          progname, filename, line, atoms[0], atoms[j]);
808 */
809                 push_bond (m, atoms[0], atoms[j]);
810               }
811         }
812       else
813         {
814           char *s1 = strdup (s);
815           for (ss = s1; *ss && *ss != '\n'; ss++)
816             ;
817           *ss = 0;
818           fprintf (stderr, "%s: %s: %d: unrecognised line: %s\n",
819                    progname, filename, line, s1);
820         }
821
822       while (*s && *s != '\n')
823         s++;
824       if (*s == '\n')
825         s++;
826       line++;
827     }
828 }
829
830
831 static int
832 parse_pdb_file (molecule *m, const char *name)
833 {
834   FILE *in;
835   int buf_size = 40960;
836   char *buf;
837   int line = 1;
838
839   in = fopen(name, "r");
840   if (!in)
841     {
842       char *buf = (char *) malloc(1024 + strlen(name));
843       sprintf(buf, "%s: error reading \"%s\"", progname, name);
844       perror(buf);
845       return -1;
846     }
847
848   buf = (char *) malloc (buf_size);
849
850   while (fgets (buf, buf_size-1, in))
851     {
852       char *s;
853       for (s = buf; *s; s++)
854         if (*s == '\r') *s = '\n';
855       parse_pdb_data (m, buf, name, line++);
856     }
857
858   free (buf);
859   fclose (in);
860
861   if (!m->natoms)
862     {
863       fprintf (stderr, "%s: file %s contains no atomic coordinates!\n",
864                progname, name);
865       return -1;
866     }
867
868   if (!m->nbonds && do_bonds)
869     {
870       fprintf (stderr, "%s: warning: file %s contains no atomic bond info.\n",
871                progname, name);
872       do_bonds = 0;
873     }
874
875   return 0;
876 }
877
878
879 typedef struct { char *atom; int count; } atom_and_count;
880
881 /* When listing the components of a molecule, the convention is to put the
882    carbon atoms first, the hydrogen atoms second, and the other atom types
883    sorted alphabetically after that (although for some molecules, the usual
884    order is different: we special-case a few of those.)
885  */
886 static int
887 cmp_atoms (const void *aa, const void *bb)
888 {
889   const atom_and_count *a = (atom_and_count *) aa;
890   const atom_and_count *b = (atom_and_count *) bb;
891   if (!a->atom) return  1;
892   if (!b->atom) return -1;
893   if (!strcmp(a->atom, "C")) return -1;
894   if (!strcmp(b->atom, "C")) return  1;
895   if (!strcmp(a->atom, "H")) return -1;
896   if (!strcmp(b->atom, "H")) return  1;
897   return strcmp (a->atom, b->atom);
898 }
899
900 static void special_case_formula (char *f);
901
902 static void
903 generate_molecule_formula (molecule *m)
904 {
905   char *buf = (char *) malloc (m->natoms * 10);
906   char *s = buf;
907   int i;
908   atom_and_count counts[200];
909   memset (counts, 0, sizeof(counts));
910   *s = 0;
911   for (i = 0; i < m->natoms; i++)
912     {
913       int j = 0;
914       char *a = (char *) m->atoms[i].label;
915       char *e;
916       while (!isalpha(*a)) a++;
917       a = strdup (a);
918       for (e = a; isalpha(*e); e++);
919       *e = 0;
920       while (counts[j].atom && !!strcmp(a, counts[j].atom))
921         j++;
922       if (counts[j].atom)
923         free (a);
924       else
925         counts[j].atom = a;
926       counts[j].count++;
927     }
928
929   i = 0;
930   while (counts[i].atom) i++;
931   qsort (counts, i, sizeof(*counts), cmp_atoms);
932
933   i = 0;
934   while (counts[i].atom)
935     {
936       strcat (s, counts[i].atom);
937       free (counts[i].atom);
938       s += strlen (s);
939       if (counts[i].count > 1)
940         sprintf (s, "[%d]", counts[i].count);  /* use [] to get subscripts */
941       s += strlen (s);
942       i++;
943     }
944
945   special_case_formula (buf);
946
947   if (!m->label) m->label = strdup("");
948   s = (char *) malloc (strlen (m->label) + strlen (buf) + 2);
949   strcpy (s, m->label);
950   strcat (s, "\n");
951   strcat (s, buf);
952   free ((char *) m->label);
953   free (buf);
954   m->label = s;
955 }
956
957 /* thanks to Rene Uittenbogaard <ruittenb@wish.nl> */
958 static void
959 special_case_formula (char *f)
960 {
961   if      (!strcmp(f, "H[2]Be"))   strcpy(f, "BeH[2]");
962   else if (!strcmp(f, "H[3]B"))    strcpy(f, "BH[3]");
963   else if (!strcmp(f, "H[3]N"))    strcpy(f, "NH[3]");
964   else if (!strcmp(f, "CHN"))      strcpy(f, "HCN");
965   else if (!strcmp(f, "CKN"))      strcpy(f, "KCN");
966   else if (!strcmp(f, "H[4]N[2]")) strcpy(f, "N[2]H[4]");
967   else if (!strcmp(f, "Cl[3]P"))   strcpy(f, "PCl[3]");
968   else if (!strcmp(f, "Cl[5]P"))   strcpy(f, "PCl[5]");
969 }
970
971
972 static void
973 insert_vertical_whitespace (char *string)
974 {
975   while (*string)
976     {
977       if ((string[0] == ',' ||
978            string[0] == ';' ||
979            string[0] == ':') &&
980           string[1] == ' ')
981         string[0] = ' ', string[1] = '\n';
982       string++;
983     }
984 }
985
986
987 /* Construct the molecule data from either: the builtins; or from
988    the (one) .pdb file specified with -molecule.
989  */
990 static void
991 load_molecules (ModeInfo *mi)
992 {
993   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
994   int wire = MI_IS_WIREFRAME(mi);
995   int i;
996
997   mc->nmolecules = 0;
998   if (molecule_str && *molecule_str && 
999       strcmp(molecule_str, "(default)"))        /* try external PDB files */
1000     {
1001       /* The -molecule option can point to a .pdb file, or to
1002          a directory of them.
1003       */
1004       struct stat st;
1005       int nfiles = 0;
1006       int list_size = 0;
1007       char **files = 0;
1008       int molecule_ctr;
1009
1010       if (!stat (molecule_str, &st) &&
1011           S_ISDIR (st.st_mode))
1012         {
1013           char buf [255];
1014           DIR *pdb_dir;
1015           struct dirent *dentry;
1016
1017           pdb_dir = opendir (molecule_str);
1018           if (! pdb_dir)
1019             {
1020               sprintf (buf, "%.100s: %.100s", progname, molecule_str);
1021               perror (buf);
1022               exit (1);
1023             }
1024
1025           if (verbose_p)
1026             fprintf (stderr, "%s: directory %s\n", progname, molecule_str);
1027
1028           nfiles = 0;
1029           list_size = 100;
1030           files = (char **) calloc (sizeof(*files), list_size);
1031
1032           while ((dentry = readdir (pdb_dir)))
1033             {
1034               int L = strlen (dentry->d_name);
1035               if (L > 4 && !strcasecmp (dentry->d_name + L - 4, ".pdb"))
1036                 {
1037                   char *fn;
1038                   if (nfiles >= list_size-1)
1039                     {
1040                       list_size = (list_size + 10) * 1.2;
1041                       files = (char **)
1042                         realloc (files, list_size * sizeof(*files));
1043                       if (!files)
1044                         {
1045                         OOM:
1046                           fprintf (stderr, "%s: out of memory (%d files)\n",
1047                                    progname, nfiles);
1048                           exit (1);
1049                         }
1050                     }
1051
1052                   fn = (char *) malloc (strlen (molecule_str) + L + 10);
1053                   if (!fn) goto OOM;
1054                   strcpy (fn, molecule_str);
1055                   if (fn[strlen(fn)-1] != '/') strcat (fn, "/");
1056                   strcat (fn, dentry->d_name);
1057                   files[nfiles++] = fn;
1058                   if (verbose_p)
1059                     fprintf (stderr, "%s: file %s\n", progname, fn);
1060                 }
1061             }
1062           closedir (pdb_dir);
1063
1064           if (nfiles == 0)
1065             fprintf (stderr, "%s: no .pdb files in directory %s\n",
1066                      progname, molecule_str);
1067         }
1068       else
1069         {
1070           files = (char **) malloc (sizeof (*files));
1071           nfiles = 1;
1072           files[0] = strdup (molecule_str);
1073           if (verbose_p)
1074             fprintf (stderr, "%s: file %s\n", progname, molecule_str);
1075         }
1076
1077       mc->nmolecules = nfiles;
1078       mc->molecules = (molecule *) calloc (sizeof (molecule), mc->nmolecules);
1079       molecule_ctr = 0;
1080       for (i = 0; i < mc->nmolecules; i++)
1081         {
1082           if (verbose_p)
1083             fprintf (stderr, "%s: reading %s\n", progname, files[i]);
1084           if (!parse_pdb_file (&mc->molecules[molecule_ctr], files[i]))
1085             {
1086               if ((wire || !do_atoms) &&
1087                   !do_labels &&
1088                   mc->molecules[molecule_ctr].nbonds == 0)
1089                 {
1090                   /* If we're not drawing atoms (e.g., wireframe mode), and
1091                      there is no bond info, then make sure labels are turned
1092                      on, or we'll be looking at a black screen... */
1093                   fprintf (stderr, "%s: %s: no bonds: turning -label on.\n",
1094                            progname, files[i]);
1095                   do_labels = 1;
1096                 }
1097               free (files[i]);
1098               files[i] = 0;
1099               molecule_ctr++;
1100             }
1101         }
1102
1103       free (files);
1104       files = 0;
1105       mc->nmolecules = molecule_ctr;
1106     }
1107
1108   if (mc->nmolecules == 0)      /* do the builtins if no files */
1109     {
1110       mc->nmolecules = countof(builtin_pdb_data);
1111       mc->molecules = (molecule *) calloc (sizeof (molecule), mc->nmolecules);
1112       for (i = 0; i < mc->nmolecules; i++)
1113         {
1114           char name[100];
1115           sprintf (name, "<builtin-%d>", i);
1116           if (verbose_p) fprintf (stderr, "%s: reading %s\n", progname, name);
1117           parse_pdb_data (&mc->molecules[i], builtin_pdb_data[i], name, 1);
1118         }
1119     }
1120
1121   for (i = 0; i < mc->nmolecules; i++)
1122     {
1123       generate_molecule_formula (&mc->molecules[i]);
1124       insert_vertical_whitespace ((char *) mc->molecules[i].label);
1125     }
1126 }
1127
1128
1129 \f
1130 /* Window management, etc
1131  */
1132 void
1133 reshape_molecule (ModeInfo *mi, int width, int height)
1134 {
1135   GLfloat h = (GLfloat) height / (GLfloat) width;
1136
1137   glViewport (0, 0, (GLint) width, (GLint) height);
1138
1139   glMatrixMode(GL_PROJECTION);
1140   glLoadIdentity();
1141   gluPerspective (30.0, 1/h, 20.0, 100.0);
1142
1143   glMatrixMode(GL_MODELVIEW);
1144   glLoadIdentity();
1145   gluLookAt( 0.0, 0.0, 30.0,
1146              0.0, 0.0, 0.0,
1147              0.0, 1.0, 0.0);
1148
1149   glClear(GL_COLOR_BUFFER_BIT);
1150 }
1151
1152
1153 static void
1154 gl_init (ModeInfo *mi)
1155 {
1156   static GLfloat pos[4] = {1.0, 0.4, 0.9, 0.0};
1157   static GLfloat amb[4] = {0.0, 0.0, 0.0, 1.0};
1158   static GLfloat dif[4] = {0.8, 0.8, 0.8, 1.0};
1159   static GLfloat spc[4] = {1.0, 1.0, 1.0, 1.0};
1160   glLightfv(GL_LIGHT0, GL_POSITION, pos);
1161   glLightfv(GL_LIGHT0, GL_AMBIENT,  amb);
1162   glLightfv(GL_LIGHT0, GL_DIFFUSE,  dif);
1163   glLightfv(GL_LIGHT0, GL_SPECULAR, spc);
1164 }
1165
1166
1167 static void
1168 startup_blurb (ModeInfo *mi)
1169 {
1170   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1171   const char *s = "Constructing molecules...";
1172   print_gl_string (mi->dpy, mc->xfont2, mc->font2_dlist,
1173                    mi->xgwa.width, mi->xgwa.height,
1174                    10, mi->xgwa.height - 10,
1175                    s);
1176   glFinish();
1177   glXSwapBuffers(MI_DISPLAY(mi), MI_WINDOW(mi));
1178 }
1179
1180 Bool
1181 molecule_handle_event (ModeInfo *mi, XEvent *event)
1182 {
1183   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1184
1185   if (event->xany.type == ButtonPress &&
1186       event->xbutton.button == Button1)
1187     {
1188       mc->button_down_p = True;
1189       gltrackball_start (mc->trackball,
1190                          event->xbutton.x, event->xbutton.y,
1191                          MI_WIDTH (mi), MI_HEIGHT (mi));
1192       return True;
1193     }
1194   else if (event->xany.type == ButtonRelease &&
1195            event->xbutton.button == Button1)
1196     {
1197       mc->button_down_p = False;
1198       return True;
1199     }
1200   else if (event->xany.type == ButtonPress &&
1201            (event->xbutton.button == Button4 ||
1202             event->xbutton.button == Button5))
1203     {
1204       gltrackball_mousewheel (mc->trackball, event->xbutton.button, 10,
1205                               !!event->xbutton.state);
1206       return True;
1207     }
1208   else if (event->xany.type == KeyPress)
1209     {
1210       KeySym keysym;
1211       char c = 0;
1212       XLookupString (&event->xkey, &c, 1, &keysym, 0);
1213
1214       if (c == ' ' || c == '\t' || c == '\r' || c == '\n')
1215         {
1216           GLfloat speed = 4.0;
1217           mc->mode = 1;
1218           mc->mode_tick = 10 * speed;
1219           return True;
1220         }
1221     }
1222   else if (event->xany.type == MotionNotify &&
1223            mc->button_down_p)
1224     {
1225       gltrackball_track (mc->trackball,
1226                          event->xmotion.x, event->xmotion.y,
1227                          MI_WIDTH (mi), MI_HEIGHT (mi));
1228       return True;
1229     }
1230
1231   return False;
1232 }
1233
1234
1235 void 
1236 init_molecule (ModeInfo *mi)
1237 {
1238   molecule_configuration *mc;
1239   int wire;
1240
1241   if (!mcs) {
1242     mcs = (molecule_configuration *)
1243       calloc (MI_NUM_SCREENS(mi), sizeof (molecule_configuration));
1244     if (!mcs) {
1245       fprintf(stderr, "%s: out of memory\n", progname);
1246       exit(1);
1247     }
1248   }
1249
1250   mc = &mcs[MI_SCREEN(mi)];
1251
1252   if ((mc->glx_context = init_GL(mi)) != NULL) {
1253     gl_init(mi);
1254     reshape_molecule (mi, MI_WIDTH(mi), MI_HEIGHT(mi));
1255   }
1256
1257   load_fonts (mi);
1258   startup_blurb (mi);
1259
1260   wire = MI_IS_WIREFRAME(mi);
1261
1262   {
1263     Bool spinx=False, spiny=False, spinz=False;
1264     double spin_speed   = 0.5;
1265     double spin_accel   = 0.3;
1266     double wander_speed = 0.01;
1267
1268     char *s = do_spin;
1269     while (*s)
1270       {
1271         if      (*s == 'x' || *s == 'X') spinx = True;
1272         else if (*s == 'y' || *s == 'Y') spiny = True;
1273         else if (*s == 'z' || *s == 'Z') spinz = True;
1274         else
1275           {
1276             fprintf (stderr,
1277          "%s: spin must contain only the characters X, Y, or Z (not \"%s\")\n",
1278                      progname, do_spin);
1279             exit (1);
1280           }
1281         s++;
1282       }
1283
1284     mc->rot = make_rotator (spinx ? spin_speed : 0,
1285                             spiny ? spin_speed : 0,
1286                             spinz ? spin_speed : 0,
1287                             spin_accel,
1288                             do_wander ? wander_speed : 0,
1289                             (spinx && spiny && spinz));
1290     mc->trackball = gltrackball_init ();
1291   }
1292
1293   orig_do_labels = do_labels;
1294   orig_do_atoms  = do_atoms;
1295   orig_do_bonds  = do_bonds;
1296   orig_do_shells = do_shells;
1297   orig_wire = MI_IS_WIREFRAME(mi);
1298
1299   mc->molecule_dlist = glGenLists(1);
1300   if (do_shells)
1301     mc->shell_dlist = glGenLists(1);
1302
1303   load_molecules (mi);
1304   mc->which = random() % mc->nmolecules;
1305
1306   mc->no_label_threshold = get_float_resource ("noLabelThreshold",
1307                                                "NoLabelThreshold");
1308   mc->wireframe_threshold = get_float_resource ("wireframeThreshold",
1309                                                 "WireframeThreshold");
1310   mc->mode = 0;
1311
1312   if (wire)
1313     do_bonds = 1;
1314 }
1315
1316
1317 /* Put the labels on the atoms.
1318    This can't be a part of the display list because of the games
1319    we play with the translation matrix.
1320  */
1321 void
1322 draw_labels (ModeInfo *mi)
1323 {
1324   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1325   int wire = MI_IS_WIREFRAME(mi);
1326   molecule *m = &mc->molecules[mc->which];
1327   int i, j;
1328
1329   if (!do_labels)
1330     return;
1331
1332   if (!wire)
1333     glDisable (GL_LIGHTING);   /* don't light fonts */
1334
1335   for (i = 0; i < m->natoms; i++)
1336     {
1337       molecule_atom *a = &m->atoms[i];
1338       GLfloat size = atom_size (a);
1339       GLfloat m[4][4];
1340
1341       glPushMatrix();
1342
1343       if (!wire)
1344         set_atom_color (mi, a, True, 1);
1345
1346       /* First, we translate the origin to the center of the atom.
1347
1348          Then we retrieve the prevailing modelview matrix (which
1349          includes any rotation, wandering, and user-trackball-rolling
1350          of the scene.
1351
1352          We set the top 3x3 cells of that matrix to be the identity
1353          matrix.  This removes all rotation from the matrix, while
1354          leaving the translation alone.  This has the effect of
1355          leaving the prevailing coordinate system perpendicular to
1356          the camera view: were we to draw a square face, it would
1357          be in the plane of the screen.
1358
1359          Now we translate by `size' toward the viewer -- so that the
1360          origin is *just in front* of the ball.
1361
1362          Then we draw the label text, allowing the depth buffer to
1363          do its work: that way, labels on atoms will be occluded
1364          properly when other atoms move in front of them.
1365
1366          This technique (of neutralizing rotation relative to the
1367          observer, after both rotations and translations have been
1368          applied) is known as "billboarding".
1369        */
1370
1371       glTranslatef(a->x, a->y, a->z);               /* get matrix */
1372       glGetFloatv (GL_MODELVIEW_MATRIX, &m[0][0]);  /* load rot. identity */
1373       m[0][0] = 1; m[1][0] = 0; m[2][0] = 0;
1374       m[0][1] = 0; m[1][1] = 1; m[2][1] = 0;
1375       m[0][2] = 0; m[1][2] = 0; m[2][2] = 1;
1376       glLoadIdentity();                             /* reset modelview */
1377       glMultMatrixf (&m[0][0]);                     /* replace with ours */
1378
1379       glTranslatef (0, 0, (size * 1.1));           /* move toward camera */
1380
1381       glRasterPos3f (0, 0, 0);                     /* draw text here */
1382
1383       /* Before drawing the string, shift the origin to center
1384          the text over the origin of the sphere. */
1385       glBitmap (0, 0, 0, 0,
1386                 -string_width (mc->xfont1, a->label) / 2,
1387                 -mc->xfont1->descent,
1388                 NULL);
1389
1390       for (j = 0; j < strlen(a->label); j++)
1391         glCallList (mc->font1_dlist + (int)(a->label[j]));
1392
1393       glPopMatrix();
1394     }
1395
1396   /* More efficient to always call glEnable() with correct values
1397      than to call glPushAttrib()/glPopAttrib(), since reading
1398      attributes from GL does a round-trip and  stalls the pipeline.
1399    */
1400   if (!wire)
1401     glEnable (GL_LIGHTING);
1402 }
1403
1404
1405 static void
1406 pick_new_molecule (ModeInfo *mi, time_t last)
1407 {
1408   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1409
1410   if (mc->nmolecules == 1)
1411     {
1412       if (last != 0) return;
1413       mc->which = 0;
1414     }
1415   else if (last == 0)
1416     {
1417       mc->which = random() % mc->nmolecules;
1418     }
1419   else
1420     {
1421       int n = mc->which;
1422       while (n == mc->which)
1423         n = random() % mc->nmolecules;
1424       mc->which = n;
1425     }
1426           
1427   if (verbose_p)
1428     {
1429       char *name = strdup (mc->molecules[mc->which].label);
1430       char *s = strpbrk (name, "\r\n");
1431       if (s) *s = 0;
1432       fprintf (stderr, "%s: drawing %s (%d)\n", progname, name, mc->which);
1433       free (name);
1434     }
1435
1436   mc->polygon_count = 0;
1437
1438   glNewList (mc->molecule_dlist, GL_COMPILE);
1439   ensure_bounding_box_visible (mi);
1440
1441   do_labels = orig_do_labels;
1442   do_atoms  = orig_do_atoms;
1443   do_bonds  = orig_do_bonds;
1444   do_shells = orig_do_shells;
1445   MI_IS_WIREFRAME(mi) = orig_wire;
1446
1447   if (mc->molecule_size > mc->no_label_threshold)
1448     do_labels = 0;
1449   if (mc->molecule_size > mc->wireframe_threshold)
1450     MI_IS_WIREFRAME(mi) = 1;
1451
1452   if (MI_IS_WIREFRAME(mi))
1453     do_bonds = 1, do_shells = 0;
1454
1455   if (!do_bonds)
1456     do_shells = 0;
1457
1458   if (! (do_bonds || do_atoms || do_labels))
1459     {
1460       /* Make sure *something* shows up! */
1461       MI_IS_WIREFRAME(mi) = 1;
1462       do_bonds = 1;
1463     }
1464
1465   build_molecule (mi, False);
1466   glEndList();
1467
1468   if (do_shells)
1469     {
1470       glNewList (mc->shell_dlist, GL_COMPILE);
1471       ensure_bounding_box_visible (mi);
1472
1473       do_labels = 0;
1474       do_atoms  = 1;
1475       do_bonds  = 0;
1476
1477       build_molecule (mi, True);
1478
1479       glEndList();
1480       do_bonds  = orig_do_bonds;
1481       do_atoms  = orig_do_atoms;
1482       do_labels = orig_do_labels;
1483     }
1484 }
1485
1486
1487 void
1488 draw_molecule (ModeInfo *mi)
1489 {
1490   static time_t last = 0;
1491   time_t now = time ((time_t *) 0);
1492   GLfloat speed = 4.0;  /* speed at which the zoom out/in happens */
1493
1494   molecule_configuration *mc = &mcs[MI_SCREEN(mi)];
1495   Display *dpy = MI_DISPLAY(mi);
1496   Window window = MI_WINDOW(mi);
1497
1498   if (!mc->glx_context)
1499     return;
1500
1501   if (last == 0)
1502     {
1503       pick_new_molecule (mi, last);
1504       last = now;
1505     }
1506   else if (mc->mode == 0)
1507     {
1508       static int tick = 0;
1509       if (tick++ > 10)
1510         {
1511           time_t now = time((time_t *) 0);
1512           if (last == 0) last = now;
1513           tick = 0;
1514
1515           if (!mc->button_down_p &&
1516               mc->nmolecules > 1 &&
1517               last + timeout <= now)
1518             {
1519               /* randomize molecules every -timeout seconds */
1520               mc->mode = 1;    /* go out */
1521               mc->mode_tick = 10 * speed;
1522               last = now;
1523             }
1524         }
1525     }
1526   else if (mc->mode == 1)   /* out */
1527     {
1528       if (--mc->mode_tick <= 0)
1529         {
1530           mc->mode_tick = 10 * speed;
1531           mc->mode = 2;  /* go in */
1532           pick_new_molecule (mi, last);
1533           last = now;
1534         }
1535     }
1536   else if (mc->mode == 2)   /* in */
1537     {
1538       if (--mc->mode_tick <= 0)
1539         mc->mode = 0;  /* normal */
1540     }
1541   else
1542     abort();
1543
1544   glPushMatrix ();
1545   glScalef(1.1, 1.1, 1.1);
1546
1547   {
1548     double x, y, z;
1549     get_position (mc->rot, &x, &y, &z, !mc->button_down_p);
1550     glTranslatef((x - 0.5) * 9,
1551                  (y - 0.5) * 9,
1552                  (z - 0.5) * 9);
1553
1554     gltrackball_rotate (mc->trackball);
1555
1556     get_rotation (mc->rot, &x, &y, &z, !mc->button_down_p);
1557     glRotatef (x * 360, 1.0, 0.0, 0.0);
1558     glRotatef (y * 360, 0.0, 1.0, 0.0);
1559     glRotatef (z * 360, 0.0, 0.0, 1.0);
1560   }
1561
1562   glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
1563
1564   if (mc->mode != 0)
1565     {
1566       GLfloat s = (mc->mode == 1
1567                    ? mc->mode_tick / (10 * speed)
1568                    : ((10 * speed) - mc->mode_tick + 1) / (10 * speed));
1569       glScalef (s, s, s);
1570     }
1571
1572   glPushMatrix();
1573   glCallList (mc->molecule_dlist);
1574
1575   if (mc->mode == 0)
1576     {
1577       molecule *m = &mc->molecules[mc->which];
1578
1579       draw_labels (mi);
1580
1581       /* This can't go in the display list, or the characters are spaced
1582          wrongly when the window is resized. */
1583       if (do_titles && m->label && *m->label)
1584         {
1585           set_atom_color (mi, 0, True, 1);
1586           print_gl_string (mi->dpy, mc->xfont2, mc->font2_dlist,
1587                            mi->xgwa.width, mi->xgwa.height,
1588                            10, mi->xgwa.height - 10,
1589                            m->label);
1590         }
1591     }
1592   glPopMatrix();
1593
1594   if (do_shells)
1595     {
1596       glColorMask (GL_FALSE, GL_FALSE, GL_FALSE, GL_FALSE);
1597       glPushMatrix();
1598       glCallList (mc->shell_dlist);
1599       glPopMatrix();
1600       glColorMask (GL_TRUE, GL_TRUE, GL_TRUE, GL_TRUE);
1601
1602       glDepthFunc (GL_EQUAL);
1603       glEnable (GL_BLEND);
1604       glBlendFunc (GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
1605       glPushMatrix();
1606       glCallList (mc->shell_dlist);
1607       glPopMatrix();
1608       glDepthFunc (GL_LESS);
1609       glDisable (GL_BLEND);
1610     }
1611
1612   glPopMatrix ();
1613
1614   mi->polygon_count = mc->polygon_count;
1615
1616   if (mi->fps_p) do_fps (mi);
1617   glFinish();
1618
1619   glXSwapBuffers(dpy, window);
1620 }
1621
1622 #endif /* USE_GL */