********> bugfix.57
Author:   Bill Ross
Date:     11/20/95

Programs: Carnal

Severity: Moderate

Problem:  Per-atom and per-residue RMS reporting is wrong / can cause crash.
          Per-atom RMS is reported as distances in the summary and as 
          distances squared in the table (fix is to use distances in
          table). Also, tables with many items (such as with rmsid%atoms) 
          can exhaust the letters used at the beginning of each line, 
          resulting in non-printing characters and duplicate letters.

Cause:    RMS: arrays and their indexes not managed correctly.
          Tables: bad design decision. New strategy: each line begins with
          L#, where '#' starts at 0 for each coordinate set analyzed.

Note:     This bugfix also removes the cosmetic problem of printing residue 
          numbers starting with 0 instead of 1 in rmsid%atoms summary output.

Fix:  
------------------------------------------------------------------------------
Make the following change to head.h:
------------------------------------------------------------------------------
*** OLD head.h
--- NEW head.h
***************
*** 175,181 ****
  			char		format[TOKLEN];
  			int		itemsline;
  			int		sofar;
! 			char		linedelim;
  			int		nfiles;
  			char		filename[TOKLEN];
  			filestruct	*files;
--- 175,181 ----
  			char		format[TOKLEN];
  			int		itemsline;
  			int		sofar;
! 			int		linedelim;
  			int		nfiles;
  			char		filename[TOKLEN];
  			filestruct	*files;
------------------------------------------------------------------------------
Make the following changes to table.c:
------------------------------------------------------------------------------
*** OLD table.c
--- NEW table.c
***************
*** 17,30 ****
  *  software.                                                           **
  *************************************************************************/
  #include "head.h"
  
  /*
!  * Table operations - 'wrap' table has a single row starting with an 'A'
!  *	that continues on the next line with a 'B' & so on; to get a
!  *	page width of columns, grep the appropriate letter from the table.
!  *	Non-wrap is intended to write a row across multiple files
   *	numbered off the table id - tabid0, tabid1...
-  *
   */
  
  /* #define DBG */
--- 17,30 ----
  *  software.                                                           **
  *************************************************************************/
  #include "head.h"
+ #include <math.h>
  
  /*
!  *  Table operations - 'wrap' table has a single row starting with 'L0'
!  *	that continues on the next line with a 'L1' & so on; to get a
!  *	page width of columns, grep the appropriate line tag from the 
!  *	table.  Non-wrap is intended to write a row across multiple files
   *	numbered off the table id - tabid0, tabid1...
   */
  
  /* #define DBG */
***************
*** 96,102 ****
  	strcpy(tabptr->idunion.tab.format, "%4e ");
  	tabptr->idunion.tab.itemsline = 10;
  	tabptr->idunion.tab.sofar = 0;
! 	tabptr->idunion.tab.linedelim = 'A';
  	tabptr->idunion.tab.print = NULL;
  	
  	if (tok[0] != ';') 
--- 96,102 ----
  	strcpy(tabptr->idunion.tab.format, "%4e ");
  	tabptr->idunion.tab.itemsline = 10;
  	tabptr->idunion.tab.sofar = 0;
! 	tabptr->idunion.tab.linedelim = 0;
  	tabptr->idunion.tab.print = NULL;
  	
  	if (tok[0] != ';') 
***************
*** 259,273 ****
  				break;
  			    case ATOMRMS:
  				printptr->p_nitems =
! 						idptr->idunion.rms.natoms;
  				printptr->p_value = (char *)
! 						&idptr->idunion.rms.atomd2;
  				break;
  			    case RESRMS:
  				printptr->p_nitems =
  						idptr->idunion.rms.nres;
  				printptr->p_value = (char *)
! 						&idptr->idunion.rms.resd;
  				break;
  			}
  			break;
--- 259,273 ----
  				break;
  			    case ATOMRMS:
  				printptr->p_nitems =
! 				  idptr->idunion.rms.grp->idunion.group.natoms;
  				printptr->p_value = (char *)
! 					idptr->idunion.rms.atomd2;
  				break;
  			    case RESRMS:
  				printptr->p_nitems =
  						idptr->idunion.rms.nres;
  				printptr->p_value = (char *)
! 						idptr->idunion.rms.resd;
  				break;
  			}
  			break;
***************
*** 340,346 ****
  {
  	printstruct	*print = tabpt->idunion.tab.print;
  	int		i, nitems;
! 	_REAL		*result;
  
  	tablereset(tabpt);
  	for (; print != NULL; print = print->next) {
--- 340,346 ----
  {
  	printstruct	*print = tabpt->idunion.tab.print;
  	int		i, nitems;
! 	_REAL		*result, value;
  
  	tablereset(tabpt);
  	for (; print != NULL; print = print->next) {
***************
*** 379,384 ****
--- 379,392 ----
  		     *  multiple-value things
  		     */
  		    case ATOMRMS:
+ 			nitems = print->p_nitems;
+ 			result = (_REAL *) print->p_value;
+ 			for (i=0; i<nitems; i++) {
+ 				value = * ((_REAL *) result+i);
+ 				value = sqrt(value);
+ 				REALtabwrite(tabpt, value);
+ 			}
+ 			break;
  		    case RESRMS:
  		    case PUCKER:
  			nitems = print->p_nitems;
***************
*** 385,391 ****
  			result = (_REAL *) print->p_value;
  			for (i=0; i<nitems; i++)
  				REALtabwrite(tabpt, 
! 					* ((_REAL *) result+i));
  			break;
  		    case GROUPCENTER:
  		    case GROUPCOM:
--- 393,399 ----
  			result = (_REAL *) print->p_value;
  			for (i=0; i<nitems; i++)
  				REALtabwrite(tabpt, 
! 					* ((_REAL *) result+i) );
  			break;
  		    case GROUPCENTER:
  		    case GROUPCOM:
***************
*** 427,433 ****
  
  	if (!(tabpt->idunion.tab.sofar % tabpt->idunion.tab.itemsline)) {
  		if (tabpt->idunion.tab.wrap) {
! 			fprintf(tabpt->idunion.tab.fptr, "\n%c ",
  						tabpt->idunion.tab.linedelim);
  			tabpt->idunion.tab.linedelim++; /* check range on setup */
  		} else {
--- 435,441 ----
  
  	if (!(tabpt->idunion.tab.sofar % tabpt->idunion.tab.itemsline)) {
  		if (tabpt->idunion.tab.wrap) {
! 			fprintf(tabpt->idunion.tab.fptr, "\nL%d ",
  						tabpt->idunion.tab.linedelim);
  			tabpt->idunion.tab.linedelim++; /* check range on setup */
  		} else {
***************
*** 460,466 ****
  idstruct	*tabpt;
  {
  	if (tabpt->idunion.tab.wrap) {
! 		tabpt->idunion.tab.linedelim = 'A';
  		tabpt->idunion.tab.sofar = 0;
  	} else {
  		tabpt->idunion.tab.curfile = tabpt->idunion.tab.files;
--- 468,474 ----
  idstruct	*tabpt;
  {
  	if (tabpt->idunion.tab.wrap) {
! 		tabpt->idunion.tab.linedelim = 0;
  		tabpt->idunion.tab.sofar = 0;
  	} else {
  		tabpt->idunion.tab.curfile = tabpt->idunion.tab.files;
------------------------------------------------------------------------------
Make the following changes to rms.c:
------------------------------------------------------------------------------
*** OLD rms.c
--- NEW rms.c
***************
*** 545,550 ****
--- 545,552 ----
  #ifdef DEBUG
  		printf("just comparing, (not fitting)..\n");
  #endif
+ 		crd2 = idptr->idunion.rms.crd2;
+ 		refcrd = idptr->idunion.rms.refcrd;
  		for (k=0; k<grpnatoms; k++) {
  			i = Xlist[k];
  			d2 = distancesq(&refcrd[i], &crd2[i]);
***************
*** 645,650 ****
--- 647,654 ----
  #ifdef DBG
  		printf("now doing per-atom (after fitting)..\n");
  #endif
+ 		crd2 = idptr->idunion.rms.crd2a;
+ 		refcrd = idptr->idunion.rms.refcrdctr;
  		for (k=0; k<grpnatoms; k++) {
  			i = Xlist[k];
  			d2 = distancesq(&refcrd[i], &crd2[i]);
------------------------------------------------------------------------------
Make the following changes to doit.c:
------------------------------------------------------------------------------
*** OLD doit.c
--- NEW doit.c
***************
*** 362,368 ****
  	idstruct	*grptr = idptr->idunion.rms.grp;
  	parmstruct	*prm = grptr->idunion.group.set->prm;
  	int		nres = idptr->idunion.rms.nres; 
! 	int		natoms = idptr->idunion.rms.natoms;
  	statstruct	*stats;
  	int		i;
  
--- 362,368 ----
  	idstruct	*grptr = idptr->idunion.rms.grp;
  	parmstruct	*prm = grptr->idunion.group.set->prm;
  	int		nres = idptr->idunion.rms.nres; 
! 	int		grpatoms = idptr->idunion.rms.grp->idunion.group.natoms;
  	statstruct	*stats;
  	int		i;
  
***************
*** 384,394 ****
  	if (idptr->idunion.rms.flag & PERATOM) {
  		printf("--RMS %s%%atoms:\n", idptr->id_name);
  		stats = idptr->idunion.rms.atomstats;
! 		for (i=0; i<natoms; i++, stats++) {
  			int	atom = atomarray[i];
  			printf("  %.4s %d  %.4s  ", 
  					&ResNames[4*AtomRes[atom]], 
! 					AtomRes[atom],
  					&AtomNames[4*atom]);
  			printstat(stats, 0.0, 0);
  		}
--- 384,394 ----
  	if (idptr->idunion.rms.flag & PERATOM) {
  		printf("--RMS %s%%atoms:\n", idptr->id_name);
  		stats = idptr->idunion.rms.atomstats;
! 		for (i=0; i<grpatoms; i++, stats++) {
  			int	atom = atomarray[i];
  			printf("  %.4s %d  %.4s  ", 
  					&ResNames[4*AtomRes[atom]], 
! 					AtomRes[atom]+1,
  					&AtomNames[4*atom]);
  			printstat(stats, 0.0, 0);
  		}
------------------------------------------------------------------------------

Temporary workarounds: None

Routines affected: rmsid%atoms and rmsid%residues, tables with
			many columns

--