//
//	kmeans3dplot.svl		3D plotting of DB fields
//
//	27-jan-2010 (rk) synchronize atom colors with $CLUSTER
//	25-jan-2010 (rk) borrowed db3dplot.svl to plot KMeans clusters
//	20-oct-2009 (jl) plot/dummy points removed on exit
//	08-oct-2009 (jl) scale graph vs. screen coordinates
//	21-aug-2007 (jd) Browse button to launch db browser
//	26-may-2006 (ms) replaced cattok with tok_cat
//	22-dec-2004 (jd) color by unique
//	10-may-2004 (jd) fixed standard deviation calculation
//	06-apr-2004 (jd) iterative monitoring loop
//	28-jul-2003 (jd) preferences panel and new thresholds
//	21-apr-2003 (jd) visibility monitor
//	10-apr-2003 (jd) selection monitor
//	08-aug-2001 (pl) fixed cov non-symmetric matrix problem
//	17-jul-2000 (pl) fixed bug when vector of numbers in cell
//	22-feb-1999 (pl) created db3dplot.svl
//
// COPYRIGHT (C) 1998-2017 CHEMICAL COMPUTING GROUP ULC  ALL RIGHTS RESERVED.
//
// PERMISSION TO USE, COPY, MODIFY AND DISTRIBUTE THIS SOFTWARE IS HEREBY
// GRANTED PROVIDED THAT: (1) UNMODIFIED OR FUNCTIONALLY EQUIVALENT CODE
// DERIVED FROM THIS SOFTWARE MUST CONTAIN THIS NOTICE; (2) ALL CODE DERIVED
// FROM THIS SOFTWARE MUST ACKNOWLEDGE THE AUTHOR(S) AND INSTITUTION(S); (3)
// THE NAMES OF THE AUTHOR(S) AND INSTITUTION(S) NOT BE USED IN ADVERTISING
// OR PUBLICITY PERTAINING TO THE DISTRIBUTION OF THE SOFTWARE WITHOUT
// SPECIFIC, WRITTEN PRIOR PERMISSION; (4) ALL CODE DERIVED FROM THIS SOFTWARE
// BE EXECUTED WITH THE MOLECULAR OPERATING ENVIRONMENT (MOE) LICENSED FROM
// CHEMICAL COMPUTING GROUP ULC
//
// CHEMICAL COMPUTING GROUP ULC DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS
// SOFTWARE, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS,
// AND IN NO EVENT SHALL CHEMICAL COMPUTING GROUP ULC BE LIABLE FOR ANY
// SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER
// RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF
// CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN
// CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
//
// Distributed by MOLSIS Inc.
//
// TO DO:
//      monitor xfields
//	remove apt in front of dbv_EntrySelect when dbv_EntrySelect vectorized
//
 
#set title	'Database k-means 3D Plot'
#set class	'MOE:interactive'
#set version	'2010.01.27'
#set main	'kmeans_Plot3D'

function FilePrompt;
function mp_SymmetricEigensystem;

const WINDOW_TITLE = 'KMeans 3D Plot';

const EXTENT = 4;
const BIN = (2*EXTENT) * (dec igen inc 40 / 40) - EXTENT;	// range

static Plot_Key;		// plot panel window key
static Pref_Axis_Data;		// pref. panel axis attr widget values
static Pref_Color_Data;		// pref. panel color widget values

static PlotAtoms;		// atoms monitored
static Monitor_ID;		// monitor task id
static Factor;			// scaling factor for plots
static Axis_e;			// e for each dimension (sqrt diagonal of cov)
static Mean;			// mean of data in each dimension
static Title;			// Titles of the axis

static MDB_file;		// name of mdb file

const TITLE_COLOR = 0x00FF00;	// default title color
const TLAB_COLOR = 0x00CC00;	// default tick label color
const GRID_COLOR = 0x00AA00;	// default grid color

const X = 1, Y = 2, Z = 3;
const e_VAL = 1, MIN_VAL = 2, MAX_VAL = 3;
const TITLE_AUTO = 1, TITLE_CUSTOM = 2;
const RANGE_SPAN = 1, RANGE_NxSD = 2, RANGE_CUSTOM = 3;
const TICK_AUTO = 1, TICK_NxSD = 2, TICK_CUSTOM = 3;
const EXCLUDE_FIELD = ['FP:*', 'PH4:*'];

// ------------------------------ GUI UTILS ----------------------------

const PLOT_DATA_WIDGET = [
	'sel_ent', 'xfield', 'yfield', 'zfield',
	'data_decorrelate', 'normalize'
    ];

const PLOT_COLOR_WIDGET = [
	'act', 'thres_val'
    ];

const PLOT_DEFAULTS = [
    normalize:		1,				// normalize the axes
    data_decorrelate:	1,				// decorrelate axes
    thres_val:		'5%, 25%, 50%, 75%, 95%'	// default thresholds
];

const WIDGET_NAME_ROOT = [ 'axis_title', 'title_pager', 'title_val',
    'axis_show_title', 'axis_show_ticks', 'axis_show_labels', 'axis_range',
    'range_pager', 'axis_N', 'axis_from', 'axis_to', 'axis_tick',
    'tick_pager', 'tick_N', 'tick_incr', 'grid_mode'];

const WIDGET_NAME_SUFFIX = ['_x', '_y', '_z'];

const PREF_AXIS_WIDGET = [
    'axis_title_x', 'axis_title_y', 'axis_title_z',
    'title_val_x', 'title_val_y', 'title_val_z',
    'axis_range_x', 'axis_range_y', 'axis_range_z',
    'axis_N_x', 'axis_N_y', 'axis_N_z',
    'axis_from_x', 'axis_from_y', 'axis_from_z',
    'axis_to_x', 'axis_to_y', 'axis_to_z',
    'axis_tick_x', 'axis_tick_y', 'axis_tick_z',
    'tick_N_x', 'tick_N_y', 'tick_N_z',
    'tick_incr_x', 'tick_incr_y', 'tick_incr_z',
    'grid_mode_x', 'grid_mode_y', 'grid_mode_z',
    'axis_show_title_x', 'axis_show_title_y', 'axis_show_title_z',
    'axis_show_ticks_x', 'axis_show_ticks_y', 'axis_show_ticks_z',
    'axis_show_labels_x', 'axis_show_labels_y', 'axis_show_labels_z'
];

const PREF_COLOR_WIDGET = [
    'grid_color_r', 'grid_color_g', 'grid_color_b',
    'tlab_color_r', 'tlab_color_g', 'tlab_color_b',
    'title_color_r', 'title_color_g', 'title_color_b'
];

const PREF_AXIS_DEFAULTS = [
    TITLE_AUTO, TITLE_AUTO, TITLE_AUTO,		// axis_title
    '', '', '',					// title_val
    RANGE_NxSD, RANGE_NxSD, RANGE_NxSD,		// axis_range
    3, 3, 3,					// axis_N
    0, 0, 0,					// axis_from
    1, 1, 1,					// axis_to
    1, 1, 1,					// axis_tick
    1, 1, 1,					// tick_N_x
    1, 1, 1,					// tick_incr
    'Backdrop', 'Backdrop', 'Backdrop',		// grid_mode
    1, 1, 1,					// axis_show_title
    1, 1, 1,					// axis_show_ticks
    1, 1, 1					// axis_show_labels
];

const MAX_SELMON_RESTART = 500;		// max number of monitor restart

// This function create a trio of RGB sliders. We append _r, _g and _b to
// parameter wname to name the sliders.

local function rgb_sliders wname
    local vbox = [
	Scale : [
	    name: tok_cat [wname,'_r'], title: 'R:',
	    range: [0,255,1], extendH: 1, len:5, onTrigger: 'return',
	    bubbleHelp: 'The red color component.'
	],
	Scale : [
	    name: tok_cat [wname,'_g'], title: 'G:',
	    range: [0,255,1], extendH: 1, len:5, onTrigger: 'return',
	    bubbleHelp: 'The green color component.'
	],
	Scale : [
	    name: tok_cat [wname,'_b'], title: 'B:',
	    range: [0,255,1], extendH: 1, len:5, onTrigger: 'return',
	    bubbleHelp: 'The blue color component.'
	]
    ];
    return vbox;
endfunction

// colors_page returns the Colors page of the Preferences panel

local function colors_page []
    return [
	extendH : 1,
	Option : [
	    title: 'Color:', name: 'color_ctrl',
	    type: 'int', onTrigger : 'return',
	    text: ['Grid', 'Tick Labels', 'Axis Title']
	],
	
	Pager : [
	    name: 'color_pager', extendH : 1,
	    Vbox : rgb_sliders 'grid_color',
	    Vbox : rgb_sliders 'tlab_color',
	    Vbox : rgb_sliders 'title_color'
	]
    ];
endfunction

// attr_page creates one Axis Attributes page of the Preferences panel using
// the widget names specified by parameter widget_name

local function attr_page widget_name
    return [
	extendV : 1,
	Hbox : [
	    extendH : 1,
	    Radio : [
		title: 'Axis Title:', name: widget_name(1), type : 'int',
		text: ['Auto', 'Custom'], onTrigger : 'return',
		bubbleHelp : [
		    'The axis is labeled using the name\n' 
		    'of the corresponding plot field.',
		    'The axis is labeled with\n'
		    'the specified text string.'
		]
	    ],
	    Pager : [
		name : widget_name(2), extendH : 1,
		Label : [],
		Text : [
		    title: 'Title:', name: widget_name(3), sensitive: 1,
		    extendH: 1, len: 20 ,
		    bubbleHelp : 
			'Title to display when the Show Title option is on.'
		]
	    ]
	],

	Hbox : [
	    Checkbox : [
		text: 'Show Title', name: widget_name(4),
		bubbleHelp :
		    'If on, the axis title is displayed.'
	    ],
	    Checkbox : [
		text: 'Show Ticks', name: widget_name(5),
		bubbleHelp :
		    'If on, axis ticks are drawn.'
	    ],
	    Checkbox : [
		text: 'Show Tick Labels', name: widget_name(6),
		bubbleHelp :
		    'If on, axis ticks are labeled with their values.'
	    ]
	],

	Hbox : [
	    extendH : 1,
	    Radio : [
		title: 'Axis Range:', name: widget_name(7), type : 'int',
		text: ['Span', 'NxSD', 'Custom'], onTrigger : 'return',
		bubbleHelp : [
		    'The axis spans the entire data range.',
		    'The axis extends N standard deviations from\n'
		    'the origin in both directions. In this case, the\n'
		    'data is shifted to have a mean of zero, so that\n'
		    'the data is centered about the origin.',
		    'The axis is drawn between the specified values.\n'
		    'The start and end points are both included.'
		]
	    ],
	    Pager : [
		name : widget_name(8), extendH : 1,
		Label : [],
		Text : [
		    title : 'N:', name: widget_name(9), len: 5,
		    type : 'real',
		    bubbleHelp : [
			'The number of standard deviations\n'
			'used to compute the axis range.'
		    ]
		],
		Hbox : [
		    Text : [
			title : 'From:', name: widget_name(10), len: 5,
			extendH : 1, type : 'real',
			bubbleHelp : [
			    'The axis start point.'
			]
		    ],
		    Text : [
			title : 'To:', name: widget_name(11), len: 5,
			extendH : 1, type : 'real',
			bubbleHelp : [
			    'The axis end point.'
			]
		    ]
		]
	    ]
	],

	Hbox : [
	    extendH : 1,
	    Radio : [
		title: 'Tick Spacing:', name: widget_name(12), type : 'int',
		text: ['Auto', 'NxSD', 'Custom'], onTrigger : 'return',
		bubbleHelp : [
		    'Automatically calculates tick spacing\n'
		    'based on the data range.',
		    'Uses N standard deviations as tick spacing unit.',
		    'Uses the specified tick spacing.'
		]
	    ],
	    Pager : [
		name : widget_name(13), extendH : 1,
		Label : [],
		Text : [
		    title : 'N:', name: widget_name(14), len: 5,
		    type : 'real',
		    bubbleHelp : [
			'The number of standard deviations\n'
			'between two adjacent ticks.'
		    ]
		],
		Text : [
		    title : 'Increment:', name: widget_name(15), len: 5,
		    type : 'real',
		    bubbleHelp : [
			'Tick spacing increment.'
		    ]
		]
	    ]
	],

	Radio : [
	    title: 'Grid Mode:', name: widget_name(16), columns: 3,
	    text: ['Backdrop', 'Stack', 'None'],
	    bubbleHelp : [
		'A plane perpendicular to the axis is placed\n'
		'at the first tick value in order to provide\n'
		'a visual cue for the data points.',
		'Planes perpendicular to the axis\n'
		'are inserted at each tick value.',
		'No grid is displayed.'
	    ]
	]
    ];
endfunction

// create_preferences_panel uses attr_page and color_page to create the
// preferences panel.

local function create_preferences_panel pref_wname
    return [
	title: 'KMeans 3D Plot Preferences', name:'pref',
	windowName: pref_wname, text: ['Apply', 'Cancel'],
	onTrigger: ['return', 'exit'],
	Option : [
	    name: 'attr_ctrl', flushLeft : 1, type: 'int',
	    text: ['X-Axis Attributes','Y-Axis Attributes',
		   'Z-Axis Attributes', 'Colors'],
	    onTrigger : 'return'
	],

	Pager : [
	    name: 'attr_pager', flushLeft : 1, extendH : 1, extendV : 1,
	    shadow : 'out', margin : 1,
	    Vbox :
		attr_page apt tok_cat [WIDGET_NAME_ROOT,WIDGET_NAME_SUFFIX(1)],
	    Vbox :
		attr_page apt tok_cat [WIDGET_NAME_ROOT,WIDGET_NAME_SUFFIX(2)],
	    Vbox :
		attr_page apt tok_cat [WIDGET_NAME_ROOT,WIDGET_NAME_SUFFIX(3)],
	    Vbox :
		colors_page []
	]
    ];
endfunction

local function create_plot_panel plot_wname
    return [
	title: WINDOW_TITLE, name:'plot', windowName: plot_wname,
	text: ['Plot', 'Browse', 'Preferences...', 'Close'],
	onTrigger: 'return',

	Text : [
	    title: 'Database:', sensitive: 0, name: 'file',
	    extendH: 1, len: 40
	],

	Checkbox : [
	    text: 'Selected Entries Only', name: 'sel_ent',
	    bubbleHelp:
		'If on, then only the selected entries in\n'
		'the Database Viewer will be plotted.'
	],

	Hbox: [
	    extendH : 1,
	    Label : [ title: 'Plot Fields:'],
	    Option : [
		title: 'X:', name: 'xfield', extendH: 1
	    ],
	    Option : [ 
		title: 'Y:', name: 'yfield', extendH: 1
	    ],
	    Option : [
		title: 'Z:', name: 'zfield', extendH: 1
	    ]
	],

	Hbox : [
	    Checkbox : [
		text: 'Decorrelate Axes', name: 'data_decorrelate',
		bubbleHelp:
		    'If on, then the three selected database fields\n'
		    'will be decorrelated prior to plotting; otherwise\n'
		    'the fields will be plotted as-is.'
	    ],
	    Checkbox : [
		text: 'Normalize Axes', name: 'normalize',
		bubbleHelp:
		    'If on, then the three selected database\n'
		    'fields will be divided by their standard\n'
		    'deviations prior to plotting.'
	    ]
	],

	Hbox : [
	    extendH : 1,
	    Option : [
		title: 'Activity:', name: 'act', onTrigger : 'return',
		bubbleHelp:
		    'The database field used to color the plotted\n'
		    'points in the 3D rendering window.'
	    ],
	    Text : [
		title: 'Threshold:', type: 'char', name: 'thres_val',
		len: 20, extendH: 1, sensitive : 0,
		shortcut : [ 
		    '5%, 25%, 50%, 75%, 95%',
		    '5%, 95%',
		    '1%...',
		    '25%...',
		    'Unique'
		],
		bubbleHelp:
		    'Thresholds used to color the data. Thresholds can\n'
		    'be specified as values or percentiles. Use the %\n'
		    'sign to indicate percentiles. Theshold values are\n'
		    'only used when an Activity field is specified.'
	    ]
	]
    ];
endfunction

// ------------------------------ GRID UTILS ----------------------------

// function tick_line is used to calculate ticks on a line. It takes for
// parameter the start and end points of the line and the tick increment.
// tick_line returns all points where a tick should be drawn. We always
// draw a tick on the start and end point

local function tick_line [val1, val2, tick]
   local tick_val = [];
   if tick > 0 then
       local minv = minE [val1, val2];
       local maxv = maxE [val1, val2];
       local range = maxv - minv;
       local tick_x = ceil [minv / tick];
       local tick_anchor = tick_x * tick;
       local ntick = floor [(maxv - tick_anchor) / tick + 1];
       tick_val = add [igen ntick, tick_x - 1] * tick;
   endif
   return uniq cat [val1, tick_val, val2];
endfunction

// function XY_Grid draw the grid for the Z-axis.

local function XY_Grid [gkey, axes, opt, grid_color]
    local i;
    local show_ticks = [opt.axis_show_ticks_x,
			opt.axis_show_ticks_y,
			opt.axis_show_ticks_z
    ];

    for i = X, Z loop
	if not show_ticks(i) then
	    axes(i) = [first axes(i), last axes(i)];
	endif
    endloop

    local len_x = length axes(X);
    local len_y = length axes(Y);
    local len_z = length axes(Z);
    local stack_len_z = len_z - 1;

	// Draw Y lines

    local npts = 2 * len_x;
    local gobj_idx = igen npts;
    local grid =  [ cat tr [axes(X), axes(X)],
		    cat rep [[first axes(Y), last axes(Y)], len_x],
		    first axes(Z)
		  ];
    GVertex [ gkey, 2, gobj_idx, grid_color, grid];

	// if in stack mode draw Y lines at every Z

    if opt.grid_mode_z == 'Stack' then
	npts = npts * stack_len_z;
	gobj_idx = igen npts;
	grid =  [   cat rep [ grid(X), stack_len_z],
		    cat rep [ grid(Y), stack_len_z],
		    cat apt rep [ dropfirst axes(Z), len_x * 2]
		];
	GVertex [ gkey, 2, gobj_idx, grid_color, grid];
    endif

	// Draw X lines

    npts = 2 * len_y;
    gobj_idx = igen npts;
    grid =  [   cat rep [[first axes(X), last axes(X)], len_y],
		cat tr [axes(Y), axes(Y)],
		first axes(Z)
	    ];
    GVertex [ gkey, 2, gobj_idx, grid_color, grid];

	// if in stack mode draw X lines at every Z

    if opt.grid_mode_z == 'Stack' then
	npts = npts * stack_len_z;
	gobj_idx = igen npts;
	grid =  [   cat rep [ grid(X), stack_len_z],
		    cat rep [ grid(Y), stack_len_z],
		    cat apt rep [ dropfirst axes(Z), len_y * 2]
		];
	GVertex [ gkey, 2, gobj_idx, grid_color, grid];
    endif
endfunction

// function XZ_Grid draw the grid for the Y-axis.

local function XZ_Grid [gkey, axes, opt, grid_color]
    local i;
    local show_ticks = [opt.axis_show_ticks_x,
			opt.axis_show_ticks_y,
			opt.axis_show_ticks_z
    ];

    for i = X, Z loop
	if not show_ticks(i) then
	    axes(i) = [first axes(i), last axes(i)];
	endif
    endloop

    local len_x = length axes(X);
    local len_y = length axes(Y);
    local len_z = length axes(Z);
    local stack_len_y = len_y - 1;

	// Draw Z lines

    local npts = 2 * len_x;
    local gobj_idx = igen npts;
    local grid =  [ cat tr [axes(X), axes(X)],
		    first axes(Y),
		    cat rep [[first axes(Z), last axes(Z)], len_x]
		  ];
    GVertex [ gkey, 2, gobj_idx, grid_color, grid];

	// if in stack mode draw Z lines at every Y

    if opt.grid_mode_y == 'Stack' then
	npts = npts * stack_len_y;
	gobj_idx = igen npts;
	grid =  [   cat rep [ grid(X), stack_len_y],
		    cat apt rep [ dropfirst axes(Y), len_x * 2],
		    cat rep [ grid(Z), stack_len_y]
		];
	GVertex [ gkey, 2, gobj_idx, grid_color, grid];
    endif

	// Draw X lines

    npts = 2 * len_z;
    gobj_idx = igen npts;
    grid =  [   cat rep [[first axes(X), last axes(X)], len_z],
		first axes(Y),
		cat tr [axes(Z), axes(Z)]
	    ];
    GVertex [ gkey, 2, gobj_idx, grid_color, grid];

	// if in stack mode draw X lines at every Y

    if opt.grid_mode_y == 'Stack' then
	npts = npts * stack_len_y;
	gobj_idx = igen npts;
	grid =  [   cat rep [ grid(X), stack_len_y],
		    cat apt rep [ dropfirst axes(Y), len_z * 2],
		    cat rep [ grid(Z), stack_len_y]
		];
	GVertex [ gkey, 2, gobj_idx, grid_color, grid];
    endif
endfunction

// function YZ_Grid draw the grid for the X-axis.

local function YZ_Grid [gkey, axes, opt, grid_color]
    local i;
    local show_ticks = [opt.axis_show_ticks_x,
			opt.axis_show_ticks_y,
			opt.axis_show_ticks_z
    ];

    for i = X, Z loop
	if not show_ticks(i) then
	    axes(i) = [first axes(i), last axes(i)];
	endif
    endloop

    local len_x = length axes(X);
    local len_y = length axes(Y);
    local len_z = length axes(Z);
    local stack_len_x = len_x - 1;

	// Draw Z lines

    local npts = 2 * len_y;
    local gobj_idx = igen npts;
    local grid =  [ first axes(X),
		    cat tr [axes(Y), axes(Y)],
		    cat rep [[first axes(Z), last axes(Z)], len_y]
		  ];
    GVertex [ gkey, 2, gobj_idx, grid_color, grid];

	// if in stack mode draw Z lines at every X

    if opt.grid_mode_x == 'Stack' then
	npts = npts * stack_len_x;
	gobj_idx = igen npts;
	grid =  [   cat apt rep [ dropfirst axes(X), len_y * 2],
		    cat rep [ grid(Y), stack_len_x],
		    cat rep [ grid(Z), stack_len_x]
		];
	GVertex [ gkey, 2, gobj_idx, grid_color, grid];
    endif

	// Draw Y lines

    npts = 2 * len_z;
    gobj_idx = igen npts;
    grid =  [   first axes(X),
		cat rep [[first axes(Y), last axes(Y)], len_z],
		cat tr [axes(Z), axes(Z)]
	    ];
    GVertex [ gkey, 2, gobj_idx, grid_color, grid];

	// if in stack mode draw Y lines at every X

    if opt.grid_mode_x == 'Stack' then
	npts = npts * stack_len_x;
	gobj_idx = igen npts;
	grid =  [   cat apt rep [ dropfirst axes(X), len_z * 2],
		    cat rep [ grid(Y), stack_len_x],
		    cat rep [ grid(Z), stack_len_x]
		];
	GVertex [ gkey, 2, gobj_idx, grid_color, grid];
    endif
endfunction

// function create_labels create from a vector of numeric tick values a
// pretty-printed vector of tokens.

local function create_labels numbers
    local i, n_numbers = length numbers, labels, a_label;
    numbers       = numbers * Factor;
    local min_val = min numbers;
    local max_len;

    if min_val < 0 then
//	labels = apt swrite ['{n:G-+}', numbers];	// line up first digit
	labels = apt swrite ['{n:G-+.3}', numbers];	// line up first digit

	    // replace the + sign by a space

	for i = 1, n_numbers loop
	    labels(i) = mput [labels(i), labels(i) == "+", " "];
	endloop
    else
//	labels = apt swrite ['{n:G-}', numbers];
	labels = apt swrite ['{n:G-.3}', numbers];
    endif

    return app token labels;
endfunction

// function tick_halfbox draw the grids for the 3 axis according to the
// specified options. The box is specified by it's two corners.

local function tick_halfbox [gkey, pt1, pt2, opt]
    local i, axes, ticks, titles;
    local text_x, text_y, text_z;
    local grid_color = bitor [ 			// get the grid color
	bitshl [ opt.grid_color_r, 16 ],	// from the panel
	bitshl [ opt.grid_color_g, 8 ],
	opt.grid_color_b
    ];
    local tlab_color = bitor [ 
	bitshl [ opt.tlab_color_r, 16 ],	// get the tick label color
	bitshl [ opt.tlab_color_g, 8 ],		// from the panel
	opt.tlab_color_b
    ];
    local title_color = bitor [ 
	bitshl [ opt.title_color_r, 16 ],	// get the title color
	bitshl [ opt.title_color_g, 8 ],	// from the panel
	opt.title_color_b
    ];
    local axis_title = [opt.axis_title_x, opt.axis_title_y, opt.axis_title_z];
    local axis_title_val = [opt.title_val_x, opt.title_val_y, opt.title_val_z];
    local axis_field = [opt.xfield, opt.yfield, opt.zfield];
    local axis_tick = [opt.axis_tick_x, opt.axis_tick_y, opt.axis_tick_z];
    local tick_N = [opt.tick_N_x, opt.tick_N_y, opt.tick_N_z];
    local tick_incr = [opt.tick_incr_x, opt.tick_incr_y, opt.tick_incr_z];
    local axis_range = [opt.axis_range_x, opt.axis_range_y, opt.axis_range_z];
    local normalize_tlabel = zero igen 3;

	// for each axis determine the tick spacing value to use

    for i = X, Z loop
	if axis_tick(i) == TICK_AUTO then
	    if axis_range(i) == RANGE_NxSD then
//		ticks(i) = Axis_e(i);
		ticks(i) = Axis_e(i) * invz Factor;
		normalize_tlabel(i) = 1;
	    else
		local range = abs [pt2(i) - pt1(i)];
		local auto_tick = range * 0.125;
//		if auto_tick < 0.5 then
//		    auto_tick = 0.5;
		if auto_tick < 0.25 * invz Factor then
		    auto_tick = 0.25 * invz Factor;
		else
		    auto_tick = round auto_tick;
		endif
		ticks(i) = auto_tick;
	    endif
	elseif axis_tick(i) == TICK_NxSD then
//	    ticks(i) = tick_N(i) * Axis_e(i);
	    ticks(i) = tick_N(i) * Axis_e(i) * invz Factor;
	    normalize_tlabel(i) = 1;
	else
//	    ticks(i) = tick_incr(i);
	    ticks(i) = tick_incr(i) * invz Factor;
	endif

	if axis_title(i) == TITLE_AUTO then
	    titles(i) = Title(i);
	else
	    titles(i) = axis_title_val(i);
	endif
    endloop

	// compute ticks, generate tick labels
 	// and calculate labels position for the X axis

    local axe_x = tick_line [pt1(X), pt2(X), ticks(X)];
    if normalize_tlabel(X) then
	text_x = create_labels (axe_x * invz Axis_e(X));
    else
	text_x = create_labels axe_x;
    endif
    local lmin_x = first axe_x - 1, lmax_x = last axe_x + 0.5;

	// compute ticks, generate tick labels
 	// and calculate labels position for the Y axis

    local axe_y = tick_line [pt1(Y), pt2(Y), ticks(Y)];
    if normalize_tlabel(Y) then
	text_y = create_labels (axe_y * invz Axis_e(Y));
    else
	text_y = create_labels axe_y;
    endif
    local lmin_y = first axe_y - 0.5, lmax_y = last axe_y + 0.5;

	// compute ticks, generate tick labels
 	// and calculate labels position for the Z axis

    local axe_z = tick_line [pt1(Z), pt2(Z), ticks(Z)];
    if normalize_tlabel(Z) then
	text_z = create_labels (axe_z * invz Axis_e(Z));
    else
	text_z = create_labels axe_z;
    endif
    local lmin_z = first axe_z - 0.5, lmax_z = last axe_z + 1;

    axes = [axe_x, axe_y, axe_z];

	// XY grid

    if opt.grid_mode_z <> 'None' then
	XY_Grid [gkey, axes, opt, grid_color];
    endif

	// XZ grid

    if opt.grid_mode_y <> 'None' then
	XZ_Grid [gkey, axes, opt, grid_color];
    endif

	// YZ grid

    if opt.grid_mode_x <> 'None' then
	YZ_Grid [gkey, axes, opt, grid_color];
    endif

	// Axis

    GVertex [ gkey, 2, igen 2, grid_color,
	[[first axes(X), inc last axes(X)], first axes(Y), first axes(Z)]];

    if opt.axis_show_title_x then
	GText [gkey, title_color, inc last axes(X) + 0.2, first axes(Y),
		    first axes(Z), titles(X)];
    endif

    GVertex [ gkey, 2, igen 2, grid_color,
	[first axes(X), [first axes(Y), inc last axes(Y)], first axes(Z)]];

    if opt.axis_show_title_y then
	GText [gkey, title_color, first axes(X), inc last axes(Y) + 0.2,
	    first axes(Z), titles(Y)];
    endif

    GVertex [ gkey, 2, igen 2, grid_color,
	[first axes(X), first axes(Y), [first axes(Z), inc last axes(Z)]]];

    if opt.axis_show_title_z then
	GText [gkey, title_color, first axes(X), first axes(Y),
	    inc last axes(Z) + 0.2, titles(Z)];
    endif

    if opt.axis_show_labels_x then
	if opt.grid_mode_z <> 'None' then
	    apt GText [gkey, tlab_color, axes(X), lmax_y, lmin_z, text_x];
	else
	    apt GText [gkey, tlab_color, axes(X), lmin_y, lmin_z, text_x];
	endif

	if opt.grid_mode_y <> 'None' then
	    apt GText [gkey, tlab_color, axes(X), lmin_y, lmax_z, text_x];
	else
	    apt GText [gkey, tlab_color, axes(X), lmin_y, lmin_z, text_x];
	endif
    endif

    if opt.axis_show_labels_y then
	if opt.grid_mode_z <> 'None' then
	    apt GText [gkey, tlab_color, lmax_x, axes(Y), lmin_z, text_y];
	else
	    apt GText [gkey, tlab_color, lmin_x, axes(Y), lmin_z, text_y];
	endif

	if opt.grid_mode_x <> 'None' then
	    apt GText [gkey, tlab_color, lmin_x, axes(Y), lmax_z, text_y];
	else
	    apt GText [gkey, tlab_color, lmin_x, axes(Y), lmin_z, text_y];
	endif
    endif

    if opt.axis_show_labels_z then
	if opt.grid_mode_x <> 'None' then
	    apt GText [gkey, tlab_color, lmin_x, lmax_y, axes(Z), text_z];
	else
	    apt GText [gkey, tlab_color, lmin_x, lmin_y, axes(Z), text_z];
	endif

	if opt.grid_mode_y <> 'None' then
	    apt GText [gkey, tlab_color, lmax_x, lmin_y, axes(Z), text_z];
	else
	    apt GText [gkey, tlab_color, lmin_x, lmin_y, axes(Z), text_z];
	endif
    endif
endfunction

// ------------------------------ PLOT UTILS ----------------------------

// sel_entry read the specified three fields values from the selected
// entries in the db. We return the selected entry keys and field values.
// We assume fields will contain three fields name tokens. We replace bad
// values by 0

local function selentry_read_col [mdb, fields]
    local nsel = dbv_nSelectedEntries mdb;
    if nsel == 0 then return []; endif

    local i = 1;
    local ekey = 0;
    local key_val = igen nsel;
    local field_val = rep [igen 3, nsel];
    local triplet;

    for i = 1, nsel loop
	while dbv_EntrySelect [mdb, ekey = db_NextEntry [mdb,ekey]] == 0 loop
	endloop
	key_val(i) = ekey;
	triplet = db_ReadFields [mdb, ekey, fields];
	(triplet | app length triplet == 0) = 0;	// remove null
	triplet = app first triplet;			// scalar only
	triplet | not finite triplet = 0;		// remove NaN, Inf
	field_val(i) = triplet;
    endloop
    field_val = tr field_val;
    return [key_val, field_val];
endfunction

// read_col reads a column out of a database replacing bad values with 0

local function read_col [mdb, field]
    local x = db_ReadColumn [mdb, field];
    (x | app length x == 0) = 0;
    x = app first x;
    x | not finite x = 0;
    return x;
endfunction

// HistogramC computes a gaussian smoothed histogram on a given set of
// values, a smoothing width and a set of bin boundaries.

local function HistogramC [x, swid, bin]
    swid = swid * sqrt 2;
    local p = app add erf (bin/swid - [x/swid]);	// gaussian smooth
    return 0.5 * (rotlpoke [p, length x] - p);
endfunction

// Distribution takes a histogram and creates a distribution in
// which each bin probability is non-zero.

local function Distribution p
    local d = add p;
    if d <= 0 then p = one p; else p = (p + inv d); endif
    return p / add p;
endfunction

// Median finds the r'th median r in [0,1], or r'th percentile, of a
// binned normalized distribution.  x are the bin values, f the density.

local function Median [x, f, r]
    local fx = f * x;
    local m = 10;

    for 50 loop
	local old_m = m;
	local a = maxE [abs (x - m), EPS_M];
	m = (add (fx/a) - (1 - 2*r)) / add (f / a);
	if abs (old_m - m) < 1e-4 then break; endif
    endloop

    return m;
endfunction

// normalize will normalize a vector of data by dividing the data by the
// standard deviation

local function normalize x
    local u = add x * invz length x;
    local s = add sqr x * invz dec length x - sqr u;
    return (x - u) * invz sqrt s;
endfunction

// DisplayData plots given data in three dimensions.  The option argument
// contains the default display flags.

local function DisplayData [ mdb, opt ]
    local axis_e = rep [1,1,1];
    local sel_ekey, pts;
    local i;

    if MOE_BATCH then return [[], axis_e]; endif

    local f = [opt.xfield, opt.yfield, opt.zfield ];
    local axis_range = [opt.axis_range_x,opt.axis_range_y,opt.axis_range_z];
    if opt.sel_ent then
	[sel_ekey, pts] = selentry_read_col [mdb, f];
    else
	pts = apt read_col [mdb, f];
    endif

    if l_length pts == 0 then return [[], axis_e]; endif

	// determine the center of the data and subtract it off
	// (this will be added back later for display, if desired)

    Mean = app add pts * invz l_length pts;	// mean is center
    pts = pts - Mean;				// center the data

	// if we are asked to decorrelate the data, we make the
	// three by three covariance matrix directly, diagonalize it,
	// and transform the data set

    local e;				// sqrt diagonal of covariance

    if opt.data_decorrelate then
	local A = [0,0,0] + [[0,0,0]];
	A(1)(1)           = add (pts(1) * pts(1));
	A(1)(2) = A(2)(1) = add (pts(1) * pts(2));
	A(1)(3) = A(3)(1) = add (pts(1) * pts(3));
	A(2)(2)           = add (pts(2) * pts(2));
	A(2)(3) = A(3)(2) = add (pts(2) * pts(3));
	A(3)(3)           = add (pts(3) * pts(3));
	A = A * invz dec l_length pts;
	[e,A] = mp_SymmetricEigensystem [A,1];
	e = maxE [0, e];
	pts = app add (A * [pts]);
    else
	e = app add sqr pts * invz dec l_length pts;
    endif

    e = sqrt e;

    if opt.normalize then
	pts = invz e * pts;
	e = one e;
    endif

	// When not using standard deviation units, do not center the
	// data points on the mean

    for i = X, Z loop
	if axis_range(i) <> RANGE_NxSD then
	    pts(i) = pts(i) + Mean(i);
	endif
    endloop

    axis_e = e;

	// assign the box corners according to the mode.  If the box
	// corners are non-null at the end of the setting we will
	// draw the bounding box of the data

    local axis_tick = [ opt.axis_tick_x, opt.axis_tick_y, opt.axis_tick_z];

	// destroy any previous system that we have created (based
	// on the chain header).  Then create the molecular system
	// that represents the data set

    local header = tok_cat [ 'KMeansPlot3D: ', db_Filename mdb ];
    oDestroy (Chains[] | cHeader Chains[] == header);
    local chain = oCreate 0;

    cSetTag[chain, ''];
    cSetName[chain, ''];
    cSetHeader[chain, header];

    local res = oCreate chain;
    rSetName [res, '*'];
    rSetType [res, 'none'];
    rSetUID [res, 1];
    rSetINS [res, " "];

    local atoms = oCreate rep[res, length pts(1)];
    aSetElement [atoms, 'Ne'];
    aSetName [atoms, 'db3d'];

	// state is related entry key

    if opt.sel_ent then
	aSetState [atoms, sel_ekey];
    else
	aSetState [atoms, db_Entries mdb];
    endif

//  aSetPos [atoms, pts]; // atom positions must be scaled to fill the view

	// Scale the graph vs. the screen resolution (in Ang)

    Factor = sqrt max [EPS_M, (add sqr axis_e / 3) * 0.5];
    aSetPos [atoms, pts * invz Factor];
    aSetInert [atoms, 1];

	// Set the static variables

    if opt.data_decorrelate then
	Title = apt tok_cat ['PCA', ['1','2','3']];
    else
	Title = f;
    endif

    PlotAtoms = atoms;
    Axis_e = axis_e;
endfunction

// read_percentiles parses a string of comma- or space-separated numbers,
// optionally terminated by a percent sign and/or followed by two or more dots
// (with or without spaces).  It returns 3 vectors:
// 1) numeric values of words correctly formatted as numbers
// 2) flags indicating which words are terminated by a percent sign
// 3) flags indicating which words are followed by two or more dots
// 4) string values of words that are not correctly formatted as numbers
//
// Example:
//	read_percentiles "3.4, 23%, abcd, 1.2.. 4.5.6,,, 7"
// returns:
//	[ [3.4,23,1.2,7], [0,1,0,0], [0,0,1,0], [ "abcd", "4.5.6" ] ]

local function read_percentiles str
    const EMARK = "\1";				// character to mark ellipsis

    if alltrue (isspace str or str == "," or str == "%") then
	return [];
    endif

	// Strip all commas

    str | str == "," = " ";			// strip commas

	// Replace all multi-dot sequences with sequences of the EMARK char

    local m = str == ".";			// mask of dots
    m = m and shiftr m;				// a dot preceded by a dot
    str | shiftl m = " ";			// replace 1st dot with space,
    str | m = EMARK;				// following dots with EMARK

	// Read string into a list of non-whitespace words (strings)

    str = first sreadb_tbl [str, 'w'];		// tokenize (into strings)

	// Prepare the mask of words that are followed by two or more dots
	// (i.e. an ellipsis).  Remove the ellipses from the list of words.

    local emask = app first str == EMARK;	// mark the ellipses
    str = str | not emask;			// remove the ellipses as words
    emask = shiftl emask | not emask;		// adjust the mask

	// Prepare the mask of words that are terminated by a percent sign.
	// Remove the percent sign from the words.

    local pmask = app last str == "%";		// mark percent signs
    str | pmask = app droplast (str | pmask);	// strip %'s

	// Read the words as real numbers.  Prepare the mask of words
	// that cause problems.

    local [val,err] = tr apt task_call [	// read as floats
	'sreadb_tbl', tr [str, 'r'], [[errmsg:'ignore']]];
    err = tok_length err;			// mark errors

	// Return the numbers and relevant masks.

    return [				// RETURN:
	val | not err,			// valid real numbers
	pmask | not err,		// percent mask of valid real numbers
	emask | not err,		// ellipsis mask of valid real numbers
	str | err			// invalid numbers (as strings)
//	val | err			// error messages for invalid numbers
    ];
endfunction

// compute the activity class of each data point.  This is
// done by various methods based on the activty data vector

local function ColorData [ mdb, opt ]
    local atoms = PlotAtoms;
    local activity = 0;
    local class, class_idx, class_seg;

    if opt.act <> '(None)' then
	activity = cat apt db_ReadFields [mdb, aState atoms, opt.act];
	(activity | app length activity == 0) = 0;	// replace null
	activity = app first activity;			// scalar only
	activity | not finite activity = 0;		// replace NaN, Inf
    else
	aSetColorBy [atoms, 'element'];
	return;
    endif

	// We use the atom charge to hold the data points activity

    aSetCharge [atoms, activity];
    local n_act = length activity;
    local [act_idx, act_seg] = sam activity;

	// sort atoms key and activity by activity

    atoms = atoms[act_idx];
    activity = activity[act_idx];
    local percentile = igen n_act * 100;
    local min_act = first activity;
    local max_act = last activity;

	// read coloring thresholds

    local thres_list, percent_m, incr_m, discard;
    if token toupper string opt.thres_val == 'UNIQUE' then
	local uniq_act = uniq activity;
	[thres_list, percent_m, incr_m, discard] =
	    [uniq activity, zero uniq_act, zero uniq_act, []];
    else
	[thres_list, percent_m, incr_m, discard] =
	    read_percentiles string opt.thres_val;
    endif

	// if the discard list is not empty then there is a syntax error

    if length discard then
	return 'Syntax error in threshold list.';
    endif

	// We need at least one threshold

    if length thres_list == 0 then
	return 'No threshold specified.';
    endif

	// No negetive percentile allowed

    if anytrue ((thres_list | percent_m) < 0) then
	return 'Negative percentiles are not allowed in threshold list.';
    endif

    if anytrue incr_m then
	local until_m = rotlpoke [percent_m, last percent_m];
	local pairs = tr [percent_m, until_m];
	pairs = pairs | incr_m;
	if anyfalse app eqL pairs then
	    return
		'Intervals can only be created between'
		' thresholds of the same type.';
	endif
    endif

    local thres_val, n_idx = 0;

	// process raw numbers threshold

    for thres_val in thres_list | not percent_m loop
	n_idx = inc n_idx;
	class_idx(n_idx) = add (activity <= thres_val);
    endloop

	// process percentile threshold

    for thres_val in (thres_list | percent_m) * n_act loop
	n_idx = inc n_idx;
	class_idx(n_idx) = add (percentile <= thres_val);
    endloop

	// compute implicits bounds for intervals

    local low_bound = 0;
    local high_bound;
    if last percent_m then
	high_bound = 100;
    else
	high_bound = max_act;
    endif

	// prev_v is all thresholds preceding the beginning of an interval
	// prev_v is used to compute the increment of the interval
	// start_v is all thresholds beginning an interval
	// end_v is all thresholds terminating an interval
	// incr_v is all intervals increment

    local prev_v = (rotrpoke [thres_list, low_bound] | incr_m);
    local start_v = (thres_list | incr_m);
    local end_v = (rotlpoke [thres_list, high_bound] | incr_m);
    local incr_v = abs (start_v - prev_v);
    local percent_flag = percent_m | incr_m;
    local i;

    for i = 1, length incr_v loop
	local start = start_v(i), end = end_v(i),  incr = incr_v(i);
	if percent_flag(i) then
	    start = start * n_act;
	    end = end * n_act;
	    incr = incr * n_act;
	    for thres_val = start + incr, end, incr loop
		n_idx = inc n_idx;
		class_idx(n_idx) = add (percentile <= thres_val);
	    endloop
	else
	    for thres_val = start + incr, end, incr loop
		n_idx = inc n_idx;
		class_idx(n_idx) = add (activity <= thres_val);
	    endloop
	endif
    endloop

	// compute classes

    class_idx = class_idx | class_idx;		// remove 0 from the idx list
    class_idx = sort uniq class_idx;
    class = zero activity;
    class[class_idx] = igen length class_idx;
    class = rotrpoke [class, 0];		// include boundary in class
    class_seg = mtoc class;
    class = smear class;
    class = inc class;

	// group atoms by class in residues

    local parent_chain = aChain first atoms;
    local old_res = uniq oParent atoms;
    local new_res = oCreate rep [parent_chain, length class_seg];
    oReparent [atoms, cat apt rep [new_res, class_seg]];
    oDestroy old_res;

	// color the atoms by class (coloring depends on method)

    if length class_seg > 1 then
	aSetScalar [atoms, dec class / dec length class_seg];
	aSetColorBy [atoms, 'scalar'];
    else
	aSetColorBy [atoms, 'element'];
    endif
endfunction

// MonitoringLoop compare the selection state of the specified atoms with
// the selection state of the related entry in a database viewer. We use
// the atom state property to bind an entry key to an atom.

local function MonitoringLoop [mdb, atoms]
    local akey, n_atoms;
    local sel_mask, vis_mask;
    local i;

    Monitor_ID = task_id [];

	// Remove points with invalid atom key

    atoms = atoms | oValid atoms;

	// Remove points with invalid entry key

    atoms = atoms[x_join [aState atoms, dbv_Entries mdb]];
    PlotAtoms = atoms;
    n_atoms = length atoms;

    if isnull atoms then	// No points left to monitor
	exit [];
    endif

	// show/select all points based on the visibility/selection
	// state of the related entry

    for i = 1, n_atoms loop
	akey = atoms(i);
	vis_mask(i) = dbv_EntryVisible [mdb, aState akey];
	aSetHidden [akey, not vis_mask(i)];
	sel_mask(i) = dbv_EntrySelected [mdb, aState akey];
	aSetSelected [akey, sel_mask(i)];
    endloop

	// show all points based on the cluster center, $CENTER,
	// state of the related entry

//    acent = db_ReadColumn [mdb, '$CENTER'];
//    apt aSetNucleusLook [
//	PlotAtoms,
//	apt select ['small-sphere', 'point', acent]
//    ];

    local avis, evis;
    local asel, esel;
    local acent, ecent;
    local acol, ecol, ocol = zero igen n_atoms;
    local prio = task_prio 3;
    loop
	for i = 1, n_atoms loop	// For each point of the plot
	    akey = atoms(i);

		// Was there any change in visibility

	    avis = not aHidden akey;
	    evis = dbv_EntryVisible [mdb, aState akey];

	    if not eqL [vis_mask(i), avis, evis] then
		vis_mask(i) = not vis_mask(i);
		aSetHidden [akey, not vis_mask(i)];
		dbv_EntrySetVisible [mdb, aState akey, vis_mask(i)];
	    endif

		// Was there any change in selection

	    asel = aSelected akey;
	    esel = dbv_EntrySelected [mdb, aState akey];

	    if not eqL [sel_mask(i), asel, esel] then
		sel_mask(i) = not sel_mask(i);
		aSetSelected [akey, sel_mask(i)];
		dbv_EntrySetSelected [mdb, aState akey, sel_mask(i)];
	    endif
	endloop

		// Was there any change in $CENTER

	acent = (aNucleusLook atoms == 'small-sphere');
	ecent = db_ReadColumn [mdb, '$CENTER'];

	if not eqL [acent, ecent] then
	    apt aSetNucleusLook [
		PlotAtoms,
		apt select ['small-sphere', 'point', ecent]
	    ];
	    ColorData [mdb, [act: '$CLUSTER', thres_val: 'Unique']];
	endif

//		// Was there any change in $CLUSTER
//
////	acol = (aColor atoms == 'small-sphere');
//	ecol = db_ReadColumn [mdb, '$CLUSTER'];
//	if not eqL [ecol, ocol] then
////	    apt aSetNucleusLook [
////		PlotAtoms,
////		apt select ['small-sphere', 'point', ecent]
////	    ];
//	    ColorData [mdb, [act: '$CLUSTER', thres_val: 'Unique']];
//	endif
//	ocol = ecol;
	sleep 0.5;
    endloop
    task_prio prio;
    exit [];
endfunction

// StartMonitor fork a new selection monitor. If the selection monitor fail
// for any reason we will try to restart it MAX_SELMON_RESTART time.

local function StartMonitor mdb
    local i = 0; 

    while i < MAX_SELMON_RESTART loop
	if dbv_Key mdb == 0 then
	   sleep 0.5;
	   continue;
	endif

	if isnull PlotAtoms then		// nothing to do
	    sleep 0.3;
	    continue;
	endif

	if second task_wfork [master:'parent',errmsg:'ignore'] == 'child' then
	    MonitoringLoop [mdb, PlotAtoms];
	endif
	i = inc i;
    endloop
    exit[];
endfunction

// ------------------------------ USER INTERFACE ----------------------------

// Draw_Ticks is the top level function used to draw the axis, titles, tick
// labels and reference grids.

local function Draw_Ticks [mdb, opt]
    if allfalse PlotAtoms then return; endif
    local i;
    local gobj = GCreateUnique tok_cat ['KMeansPlot3D: ', db_Filename mdb];
    local axis_range = [opt.axis_range_x,opt.axis_range_y,opt.axis_range_z];
    local old_axis_range = Pref_Axis_Data[
	indexof [['axis_range_x', 'axis_range_y', 'axis_range_z'],
	    PREF_AXIS_WIDGET
	]
    ];
    local axis_Nx = abs [opt.axis_N_x, opt.axis_N_y, opt.axis_N_z];
    local axis_boundary = [
	[opt.axis_from_x, opt.axis_to_x],
	[opt.axis_from_y, opt.axis_to_y],
	[opt.axis_from_z, opt.axis_to_z]
    ];
    local corner1, corner2;

    local pos = aPos (PlotAtoms | oValid PlotAtoms);
    if not l_length pos then return; endif

    for i = X, Z loop
	if axis_range(i) <> old_axis_range(i) then
	    pos(i) = pos(i) * Factor;
	    if axis_range(i) == RANGE_NxSD then
		pos(i) = pos(i) - Mean(i);
	    elseif old_axis_range(i) == RANGE_NxSD then
		pos(i) = pos(i) + Mean(i);
	    endif
	    pos(i) = pos(i) * invz Factor;
	    aSetPos [PlotAtoms | oValid PlotAtoms, pos];
	endif
    endloop

	// for each axis determine min and max values.

    for i = X, Z loop
	if axis_range(i) == RANGE_NxSD then
	    local n = axis_Nx(i);
//	    corner1(i) = neg n * Axis_e(i);
//	    corner2(i) =  n * Axis_e(i);
	    corner1(i) = neg n * Axis_e(i) * invz Factor;
	    corner2(i) =  n * Axis_e(i) * invz Factor;
	elseif axis_range(i) == RANGE_SPAN then
	    corner1(i) = minE pos(i);
	    corner2(i) = maxE pos(i);
	elseif axis_range(i) == RANGE_CUSTOM then
//	    corner1(i) = minE axis_boundary(i);
//	    corner2(i) = maxE axis_boundary(i);
	    corner1(i) = minE axis_boundary(i) * invz Factor;
	    corner2(i) = maxE axis_boundary(i) * invz Factor;
	endif
    endloop

    tick_halfbox [gobj, corner1, corner2, opt];
endfunction

// Preference panel window loop

local function Pref_WindowLoop [pref_wname, mdb, fields]

    if WindowShow [pref_wname, 1] then return; endif

    local i, opt;
    local axis_range = ['axis_range_x', 'axis_range_y', 'axis_range_z'];
    local range_pager = ['range_pager_x','range_pager_y','range_pager_z'];
    local axis_tick = ['axis_tick_x', 'axis_tick_y', 'axis_tick_z'];
    local tick_pager = ['tick_pager_x', 'tick_pager_y', 'tick_pager_z'];
    local axis_title = ['axis_title_x', 'axis_title_y', 'axis_title_z'];
    local title_pager = ['title_pager_x', 'title_pager_y', 'title_pager_z'];
    local color_slider = ['grid_color*', 'tlab_color*', 'title_color*'];

    local pref_key = WindowCreate create_preferences_panel [pref_wname];
    WindowSetData [pref_key, tag [PREF_AXIS_WIDGET, Pref_Axis_Data]];
    WindowSetData [pref_key, tag [PREF_COLOR_WIDGET, Pref_Color_Data]];
    WindowSetAttr [pref_key, tag [range_pager, nest [page:2]]];
    WindowShow [pref_key, 1];

    loop
	local [v,trig] = WindowWait pref_key;

	if trig === 'pref' then
	    if v.pref == 'Apply' then
		opt = cat [v, WindowValues Plot_Key];
		Draw_Ticks [mdb, opt];

		Pref_Axis_Data = last untag
		    WindowGetData [pref_key, PREF_AXIS_WIDGET];
		Pref_Color_Data = last untag
		    WindowGetData[pref_key, PREF_COLOR_WIDGET];
	    elseif v.pref == 'Cancel' then
		WindowShow [pref_key, 0];
		WindowDestroy pref_key;
		break;
	    endif
	endif

	    // XYZ-axis, colors pulldown

	if trig === 'attr_ctrl' then
	    WindowSetAttr [pref_key, ['attr_pager' : [page : v.(trig)]]];
	    continue;
	endif

	    // colors page pulldown (grid, tick labels, axis title)

	if trig === 'color_ctrl' then
	    WindowSetAttr [pref_key, ['color_pager' : [page : v.(trig)]]];
	    continue;
	endif

	    // range radio

	i = indexof [trig, axis_range];
	if anytrue i then
	    WindowSetAttr [pref_key,
		tag [range_pager(i), nest [page : v.(trig)]]];
	    continue;
	endif

	    // tick spacing radio

	i = indexof [trig, axis_tick];
	if anytrue i then
	    WindowSetAttr [pref_key,
		tag [tick_pager(i), nest [page : v.(trig)]]];
	    continue;
	endif

	    // title radio

	i = indexof [trig, axis_title];
	if anytrue i then
	    WindowSetAttr [pref_key,
		tag [title_pager(i), nest [page : v.(trig)]]];
	    continue;
	endif

	    // Color sliders are echoed immediately to main window
	    // to facilitate color selection

	if m_findmatch [color_slider, trig] then
	    Pref_Color_Data = last untag
		WindowGetData[pref_key, PREF_COLOR_WIDGET];

	    opt = cat [ WindowValues Plot_Key,
		      tag [PREF_AXIS_WIDGET, Pref_Axis_Data],
		      tag [PREF_COLOR_WIDGET, Pref_Color_Data]];
	    Draw_Ticks [mdb, opt];
	endif
    endloop
endfunction

// Plot panel window loop

local function _kmeans_Plot3D [mdb_file, opt]
    if MOE_BATCH then return; endif

	// if we are given no file, then prompt the user to
	// select a database (note the the database need not be open
	// in a database viewer)

    if not MOE_BATCH then
	if mdb_file === [] then
	    mdb_file = task_getenv 'DB_VIEW';
	endif
	if mdb_file === [] then
	    mdb_file = FilePrompt [
		'3D Plot: Select a Database','open','*.mdb'
	    ];
	    if mdb_file === [] then return; endif
	endif
    endif

    MDB_file = mdb_file = db_Filename mdb_file;

    local plot_wname = token swrite ['kmeansPlot3D:{}', mdb_file];
    if WindowShow [plot_wname, 1] then return; endif

    local pref_wname = token swrite ['kmeansPlot3DPref:{}', mdb_file];
    local sel_wname = token swrite ['kmeansPlot3DSel:{}', mdb_file];

    Pref_Axis_Data = PREF_AXIS_DEFAULTS;
    Pref_Color_Data = [
	    bitshr [bitand [GRID_COLOR, 0xFF0000], 16],		// grid rgb
	    bitshr [bitand [GRID_COLOR, 0x00FF00], 8],
	    bitand [GRID_COLOR, 0x0000FF],
	    bitshr [bitand [TLAB_COLOR, 0xFF0000], 16],		// label rgb
	    bitshr [bitand [TLAB_COLOR, 0x00FF00], 8],
	    bitand [TLAB_COLOR, 0x0000FF],
	    bitshr [bitand [TITLE_COLOR, 0xFF0000], 16],	// title rgb
	    bitshr [bitand [TITLE_COLOR, 0x00FF00], 8],
	    bitand [TITLE_COLOR, 0x0000FF]
    ];

    local mdb = db_Open mdb_file;
    local ent = db_Entries mdb;
    local idx = igen length ent;

    Plot_Key = WindowCreate create_plot_panel plot_wname;

	// SetDBFields updates the panel items relating to the database
	// numeric fields

    function SetDBFields [];
	local option_name = ['xfield', 'yfield', 'zfield'];
	local fields = db_NumericFields mdb;
	if length fields == 0 then
	    exit tok_cat [db_Filename mdb, ':\nNo numeric fields to plot.'];
	endif
	local v;

	fields = fields | not m_findmatch [ EXCLUDE_FIELD, fields ];
	WindowSetAttr [Plot_Key, [
	    act:    [ text: prepend [ fields, '(None)']],
	    xfield: [ text: fields ],
	    yfield: [ text: fields ],
	    zfield: [ text: fields ]
	]];

	    // start pulldowns with selected fields if any

	if dbv_Key mdb then
	    local field_x = indexof [fields, first db_Fields mdb];
	    fields = get [fields, x_sort
		get[not dbv_FieldSelection mdb, field_x]];
	endif

	v = tag [ option_name, fields [minE [ length fields, [1,2,3]]]];
	v = tagpoke [v, 'act', '(None)'];
	v = tagcat [[
	    xfield: opt.xfield, yfield: opt.yfield, zfield: opt.zfield,
	    act: opt.act, thres_val: opt.thres_val
	], v];
	WindowSetData [Plot_Key, v];
	return fields;
    endfunction

	// set the initial data for the panel and show the panel
	// to the user.

    WindowSetData [Plot_Key, [ file: db_Filename mdb ]];
    WindowSetData [Plot_Key, PLOT_DEFAULTS];
    local fields = SetDBFields [];
    WindowShow Plot_Key;

    local thres_prev = last untag WindowGetData [Plot_Key, 'thres_val'];

	// Record plot and pref panels data

    local plot_data_value = second untag
	WindowGetData [Plot_Key, PLOT_DATA_WIDGET];
    local plot_color_value = second untag
	WindowGetData [Plot_Key, PLOT_COLOR_WIDGET];

	// Start monitor

    if second task_fork [master:'parent', idle:1] == 'child' then
	StartMonitor mdb;
	exit [];
    endif

	opt = tagcat [
	    opt,
	    tag [PREF_AXIS_WIDGET, Pref_Axis_Data],
	    tag [PREF_COLOR_WIDGET, Pref_Color_Data],
	    PLOT_DEFAULTS
	];

	DisplayData [mdb, opt];
	Draw_Ticks [mdb, opt];
	aSetNucleusLook [
	    PlotAtoms[idx | db_ReadColumn [mdb, '$CENTER']],
	    'small-sphere'
	];
	local msg = ColorData [mdb, opt];
	if length msg then
	    Warning msg;
	endif

	// Start plot window event loop

    loop
	local [v,trig] = WindowWait Plot_Key;

	if trig === 'plot' and v.plot === 'Close' then break; endif

	    // Start pref window event loop when apply is pressed

	if trig === 'plot' and v.plot === 'Preferences...' then
	    if second task_fork [ master : 'parent' ] == 'child' then
		Pref_WindowLoop [ pref_wname, mdb, fields ];
		exit [];
	    endif
	    continue;
	endif

	if trig === 'plot' and v.plot === 'Browse' then
	    local [result, code] =
		task_fork [master: 'parent', errmsg: 'ignore'];

		// parent

	    if code == '' then continue; endif

		// child

	    run ['dbbrowse.svl', [mdb, [use_sel: 1, subject: 'first_mol2D']]];
	    exit [];
	endif

	    // desactivate threshold list when no activity field selected

	if trig == 'act' then
	    if v.act == '(None)' then
		WindowSetAttr [ Plot_Key, [thres_val : [ sensitive : 0]]];
	    else
		WindowSetAttr [ Plot_Key, [thres_val : [ sensitive : 1]]];
	    endif
	    continue;
	endif

	    // fork a task to handle plot button press

	if second task_fork [master:'parent'] <> 'child' then
	    continue;
	endif

	WindowSetAttr [Plot_Key, [ plot : [ sensitive : [0,0,0,1]]]];

	opt = tagcat [
	    v,
	    tag [PREF_AXIS_WIDGET, Pref_Axis_Data],
	    tag [PREF_COLOR_WIDGET, Pref_Color_Data],
	    opt
	];

	if opt.sel_ent and not dbv_nSelectedEntries mdb then
	    Warning 'There are no selected entries.';
	    WindowSetAttr [Plot_Key, [ plot : 
		[sensitive : [1,1,1,1]]]];
	    exit [];
	endif

	    // stop monitor

	if anytrue Monitor_ID then
	    local prio = task_prio 0;
	    task_kill Monitor_ID;
	    PlotAtoms = [];
	    task_prio prio;
	endif

	DisplayData [mdb, opt];
	Draw_Ticks [mdb, opt];
	aSetNucleusLook [
	    PlotAtoms[idx | db_ReadColumn [mdb, '$CENTER']],
	    'small-sphere'
	];
//	local msg = ColorData [mdb, opt];
	msg = ColorData [mdb, opt];
	if length msg then
	    Warning msg;
	endif

	WindowSetAttr [Plot_Key, [ plot : 
	    [sensitive : [1,1,1,1]]]];
	exit[];
    endloop

    WindowDestroy Plot_Key;
    db_Close mdb;
endfunction

function kmeans_Plot3D [mdb_file, opt]
    write 'kmeans_Plot3D '; print [mdb_file, opt];
    if second task_fork [statics: 'new'] <> 'child' then return; endif
    MDB_file = [];	// reset in _kmeans_Plot3D [];

    if second task_wfork [idle: 1] <> 'child' then
	if not isnull MDB_file then
	    GDestroy tok_cat ['KMeansPlot3D: ', MDB_file];

	    local header = tok_cat ['KMeansPlot3D: ', MDB_file];
	    oDestroy (Chains [] | cHeader Chains [] == header);
	endif

	exit [];
    endif

    task_settitle [-1, WINDOW_TITLE];
    _kmeans_Plot3D [mdb_file, opt];

	// keep plot in MOE window
//    exit [];
endfunction

#eof