CCL Home Page
Up Directory CCL mmio_convert
/* Copyright 1995, Columbia University, all rights reserved.
 * Permission is granted to utilize and disseminate this code or
 *  document without charge, provided that (1) this copyright notice is 
 *  not removed, and (2) all changes made by other than members of the 
 *  MacroModel Development Group at Columbia University are, if further 
 *  disseminated, (a) noted as such; for example, by means of source-code 
 *  comment lines, and (b) communicated back to the author for possible 
 *  inclusion in subsequent versions. */

/****************************************************************************
 * pss. 2/95.  Program to read an mmod file, compressed or otherwise, and
 *  attempt to output as compressed a file as is possible without
 *  reordering atoms or bond tables.  Also can uncompress.
 *
 * The term "comp" is generally used in variable names and comments to
 *  refer to data stored in a compressed CT. The term "conn" is used
 *  to refer to the data that is stored in a full CT, but not a compressed
 *  CT.
 *
 * $RCSfile: mmio_convert.c,v $
 * $Revision: 1.10 $
 * $Date: 1996/01/24 17:51:58 $
 ***************************************************************************/

#if defined( GETOPT_STDLIB )
	#include 
#elif defined( GETOPT_UNISTD )
	#include 
#elif defined( GETOPT_GETOPT )
	#include 
#endif

#include 
#include 
#include 
#if defined( DBMALLOC )
	#include 
#endif
#include "mmio.h"

/* variable #definitions: */
#define FALSE 0
#define TRUE 1

/*    ...used in declarations of arrays to be malloc()ed: */
#define AR *

/* macro #definitions: */
/*    ...automagic exit with appropriate message if error is encountered: */
#define err_check( status ) \
	if( status == MMIO_ERR )  { \
		fprintf( stderr, \
		 "!!! mmio_uncompress: exiting, file %s, line %d\n", \
		 __FILE__, __LINE__ ); \
		return -1; \
	}

/* structure declarations: */
/*    ...info readable from full CT but not from compressed CT: */
struct ct_conntag {
	struct conn_atomtag {
		/* atom data that don't change with coordinates: */
		int itype;	/* mmod atom type */
		int nbond;	/* number of bonded neighbors (n's) */
		int bond_atom[ MMIO_MAXBOND ];	/* atom numbers of bonded n's */
		int bond_order[ MMIO_MAXBOND ];	/* bond orders to bonded n's */
		float charge1;			/* 1st charge field in .dat */
		float charge2;			/* 1st charge field in .dat */
		char chain;			/* pdb chain designator */
		int resnum;			/* pdb res number */
		char resname1;			/* 1-char res name */
		char resname4[ MMIO_S_STRLEN + 1 ]; /* 4-char pdb res name */
		char pdbname[ MMIO_S_STRLEN + 1 ];  /* 4-char pdb atom name */
	} AR atom;
	int natom;		  /* no. of atoms in this molecule */
};

/*    ...all info readable from compressed CT: */
struct ct_comptag {
	struct comp_atomtag {
		/* data per atom: */
		float xyz[ 3 ];	/* coords */
		int color;	/* color */
		int mmod_iatom;	/* MacroModel atom number -- index origin 1 */
	} AR atom;
	int natom;		/* number of atoms stored in this struct */
	char title[ MMIO_L_STRLEN + 1 ];	/* title from input file */
};

/* definitions of static variables: */
static int verbose_flag = FALSE;
struct ct_conntag global_conn[ 2 ];
struct ct_comptag global_comp[ 2 ];

/* prototypes: */
int main( int narg, char *arg[] );
static void verbose( char *fmt, ... );
static void usage( void );
static int diff_conn( struct ct_conntag *conn1, struct ct_conntag *conn2 );
static int get_ct( int iread, int natom, char *title,
 struct ct_conntag *newconn, struct ct_comptag *newcomp );
static int put_ct( int iwrite, struct ct_conntag *conn,
 struct ct_comptag *oldcomp, struct ct_comptag *newcomp );
static int filldiff( int AR diff, struct ct_comptag *oldcomp,
 struct ct_comptag *newcomp );
/************************************************************************/
int main( int narg, char *arg[] )
{
	int status, ct_type;
	struct ct_conntag *oldconn, *newconn, *tmpconn;
	struct ct_comptag *oldcomp, *newcomp, *tmpcomp;
	int ict= -1, natom;
	char title[ MMIO_L_STRLEN ];
	char readfname[200], writefname[200];
	int ct_type_requested = MMIO_COMPRESSED;
	int opt;
	extern int optind;
	int iread, iwrite;

	/* initializations: */
	global_conn[0].atom = global_conn[1].atom = NULL;
	global_comp[0].atom = global_comp[1].atom = NULL;
	global_conn[0].natom = global_conn[1].natom = global_comp[0].natom
	 = global_comp[1].natom = 0;
	oldconn = global_conn + 0;
	newconn = global_conn + 1;
	oldcomp = global_comp + 0;
	newcomp = global_comp + 1;

	/* parse cmdline options: */
	while( (opt=getopt(narg,arg,"cfv")) != -1 ) {
		switch( opt ) {
			case 'c':
				/* write compressed CT's where possible: */
				ct_type_requested = MMIO_COMPRESSED;
				break;
			case 'f':
				/* always write full CT's: */
				ct_type_requested = MMIO_FULL;
				break;
			case 'v':
				/* put play-by-play on stderr: */
				verbose_flag = TRUE;
				break;
			default:
				fprintf( stderr, "!!! mmio_convert main: found"
				 " illegal flag '%c'\n", opt );
				usage();
				return -1;
		}
	}

	/* see how many extra args there are, after optional flags: */
	switch( narg-optind ) {
		case 0:
			/* stdin for input, stdout for output: */
			strcpy( readfname, "-" );
			strcpy( writefname, "-" );
			break;
		case 1:
			/* arg for input fname, stdout for output: */
			strcpy( readfname, arg[optind] );
			strcpy( writefname, "-" );
			break;
		case 2:
			/* args for input & output fnames: */
			strcpy( readfname, arg[optind++] );
			strcpy( writefname, arg[optind++] );
			break;
		default:
			/* huh? */
			fprintf( stderr, "!!! mmio_convert main: found %d"
			 " args after flags\n", narg-optind );
			usage();
			return -1;
	}

	/* tell the user what he asked for: */
	verbose( "mmio_convert main: output format= %s\n",
	 mmio_return_code(ct_type_requested) );
	verbose( "mmio_convert main: readfname= %s\n",
	 strcmp("-",readfname)==0?"":readfname );
	verbose( "mmio_convert main: writefname= %s\n",
	 strcmp("-",writefname)==0?"":writefname );

	/* tell mmio library to direct its errmsgs (if any) to stderr: */
	mmio_errfile( stderr );
	
	/* open file for reading: */
	status = mmio_open( &iread, readfname, MMIO_READ );
	verbose( "mmio_convert main: mmio_open() returns %s for readfile\n",
	 mmio_return_code(status) );
	err_check( status );
	
	/* open file for writing: */
	status = mmio_open( &iwrite, writefname, MMIO_WRITE );
	verbose( "mmio_convert main: mmio_open() returns %s for writefile\n",
	 mmio_return_code(status) );
	err_check( status );

	/* get CT's from input file until EOF is detected: */
	while( (ct_type=mmio_get_ct(iread,ct_type_requested,&natom,title))
	 != MMIO_EOF ) {
		ict += 1;
		verbose( "mmio_convert main: mmio_get_ct() returns %s "
		 "for ict %d; natom= %d, title=\n'%s'\n",
		 mmio_return_code(ct_type), ict, natom, title );
		err_check( ct_type );

		if( ct_type == MMIO_FULL ) {
			/* read full CT into "new" data structures: */
			status = get_ct( iread, natom, title,
			 newconn, newcomp );
			verbose( "mmio_convert main: get_ct() returns %s\n",
			 mmio_return_code(status) );
			err_check( status );
			if( ict==0 || ct_type_requested==MMIO_FULL
			 || diff_conn(newconn,oldconn)) {
				/* this is either the first CT, or the
				 *  connectivity differs from the previous
				 *  CT, or we have requested full CT's;
				 *  write out the full CT we just read: */
				status = put_ct(
				 iwrite, newconn, newcomp, NULL );
				verbose( "mmio_convert main: put_ct() returns"
				 " %s\n", mmio_return_code(status) );
				err_check( status );
				
				/* now make the CT we just read the "old"
				 *  CT, and prepare to read the next one
				 *  into the "new" storage: */
				tmpconn = oldconn;
				oldconn = newconn;
				newconn = tmpconn;
				tmpcomp = oldcomp;
				oldcomp = newcomp;
				newcomp = tmpcomp;
				continue;
			}
			/* If the conn part of the new CT is the same
			 *  as that of the old one, and we have requested
			 *  compressed output, then at this moment the conn 
			 *  is in oldconn and the comp is in newcomp. */
		} else if( ct_type == MMIO_COMPRESSED ) {
			/* read compressed data into newcomp;  the conn data
			 *  remains untouched in oldconn: */
			status = get_ct( iread, natom, title, NULL, newcomp );
			verbose( "mmio_convert main: get_ct() returns %s\n",
			 mmio_return_code(status) );
			err_check( status );
		} else {
			/* unexpected ct_type */
			fprintf( stderr, "!!! mmio_uncompress: mmio_get_ct()"
			 " returns unexpected status %s\n",
			 mmio_return_code(status) );
			err_check( MMIO_ERR );
		}

		/* write out a new CT, taking comp info from newcomp
		 *  where it overrides the oldcomp info: */
		status = put_ct( iwrite, NULL, oldcomp, newcomp );
		verbose( "mmio_convert main: put_ct() returns"
		 " %s\n", mmio_return_code(status) );
		err_check( status );
	}

	/* this output covers MMIO_EOF returns from mmio_get_ct(): */
	verbose( "mmio_convert main: mmio_get_ct() returns %s "
	 "for ict %d\n", mmio_return_code(ct_type), ict );

	ict += 1;

	status = mmio_cleanup();
	verbose( "mmio_convert main: mmio_cleanup() returns"
	 " %s\n", mmio_return_code(status) );
	err_check( status );

	/* announce the number of CT's we processed: */
	verbose( "mmio_convert main: %d ct's found in input file\n", ict );

	return ict;
}
/************************************************************************/
static void verbose( char *fmt, ... )
 /* if verbose_flag is TRUE, print msg to stderr; otherwise, do nothing: */
{
	char buf[ 200 ];
	va_list args;

	if( verbose_flag == FALSE ) {
		return;
	}

	va_start( args, fmt );
	( void )vsprintf( buf, fmt, args );
	va_end( args );

	fprintf( stderr, "%s", buf );
}
/************************************************************************/
static void usage( void )
 /* simply print usage message: */
{
	fprintf( stderr, "Usage: mmio_convert [-cfv] [infile [outfile]]\n" );
}
/************************************************************************/
static int diff_conn( struct ct_conntag *conn1, struct ct_conntag *conn2 )
 /* return TRUE if conn1 and conn2 are different; FALSE otherwise: */
{
	int iatom, natom, ibond, nbond;
	struct conn_atomtag *conn1_atom, *conn2_atom;

	/* the CT's are different if they have different numbers of atoms: */
	if( (natom=conn1->natom) != conn2->natom ) {
		return TRUE;
	}


	/* compare relevant data for each atom: */
	for( iatom=0; iatomatom + iatom;
		conn2_atom = conn2->atom + iatom;
		/* the CT's are different if any atom has different numbers of
		 *  bonds in the two CT's: */
		if( (nbond=conn1_atom->nbond) != conn2_atom->nbond ) {
			return TRUE;
		}
		for( ibond=0; ibondbond_atom[ibond]
			 != conn2_atom->bond_atom[ibond]
			 || conn1_atom->bond_order[ibond]
			 != conn2_atom->bond_order[ibond] ) {
				return TRUE;
			}
		}
	}
	return FALSE;
}
/************************************************************************/
static int get_ct( int iread, int natom, char *title,
 struct ct_conntag *conn, struct ct_comptag *comp )
 /* call mmio routines to read conn and comp data into structures conn
  *  and comp.
  * if conn is NULL, then we're expecting compressed info, and we fill
  *  only conn: */
{
	int iatom, status;
	struct conn_atomtag *connatom;
	struct comp_atomtag *compatom;

	/* fill statically declared comp variables: */
	comp->natom = natom;
	strcpy( comp->title, title );

	/* malloc space in comp for atom data: */
	if( comp->atom == NULL ) {
		comp->atom = ( struct comp_atomtag * )malloc(
		 natom*sizeof(struct comp_atomtag) );
	} else {
		comp->atom = ( struct comp_atomtag * )realloc( comp->atom,
		 natom*sizeof(struct comp_atomtag) );
	}
	if( comp->atom == NULL ) {
		fprintf( stderr,
		 "!!! get_ct: malloc() fails for comp, natom= %d\n", natom );
		return MMIO_ERR;
	}

	/* fill comp; if conn is needed as well, malloc storage in it for
	 *  atoms, and fill it, too: */
	if( conn == NULL ) {
		/* no conn data needed (compressed CT being passed to us): */
		/*    ...get data atom by atom: */
		for( iatom=0; iatomatom + iatom;
			status = mmio_get_atom( iread,
			 &compatom->mmod_iatom, NULL, NULL, NULL, NULL,
			 compatom->xyz, NULL, NULL, NULL,
			 &compatom->color, NULL, NULL, NULL, NULL );
			verbose( "get_ct: mmio_get_atom() returns %s,"
			 " iatom= %d\n", mmio_return_code(status), iatom );
			if( status == MMIO_ERR ) break;
		}
	} else {
		/* conn data needed (full CT being passed to us): */
		/*    ...fill statically declared conn variables: */
		conn->natom = natom;

		/*    ...malloc space in conn for atom data: */
		if( conn->atom == NULL ) {
			conn->atom = ( struct conn_atomtag * )malloc(
			 natom*sizeof(struct conn_atomtag) );
		} else {
			conn->atom = ( struct conn_atomtag * )realloc(
			 conn->atom, natom*sizeof(struct conn_atomtag) );
		}
		if( conn->atom == NULL ) {
			fprintf( stderr, "!!! get_ct: malloc() fails for conn,"
			 " natom= %d\n", natom );
			return MMIO_ERR;
		}

		/*    ...get data atom by atom: */
		for( iatom=0; iatomatom + iatom;
			connatom = conn->atom + iatom;
			status = mmio_get_atom( iread,
			 &compatom->mmod_iatom, &conn->atom[iatom].itype, 
			 &connatom->nbond, connatom->bond_atom,
			 connatom->bond_order, compatom->xyz,
			 &connatom->charge1, &connatom->charge2,
			 &connatom->chain, &compatom->color,
			 &connatom->resnum, &connatom->resname1,
			 connatom->resname4, connatom->pdbname );
			verbose( "get_ct: mmio_get_atom() returns %s,"
			 " iatom= %d\n",
			 mmio_return_code(status), iatom );
			if( status == MMIO_ERR ) break;
		}
	}
	return status;
}
/************************************************************************/
static int put_ct( int iwrite, struct ct_conntag *conn,
 struct ct_comptag *oldcomp, struct ct_comptag *newcomp )
 /* call mmio routines to write out a CT.  
  * If conn is NULL, then we will be writing a compressed CT;  in this event,
  *  we write lines only for those atoms in newcomp which have coordinates 
  *  or colors differing from those of the corresponding atoms in oldcomp.
  * If conn is not NULL, then we'll write a full CT, using connectivity data
  *  from conn and coords and colors from oldcomp. */
{
	int status, natom, iatom, idiff;
	struct conn_atomtag *conn_atom;
	struct comp_atomtag *comp_atom;
	static int AR diff=NULL, ndiff=0;
	if( conn == NULL ) {
		/* write compressed CT, only for atoms where newcomp
		 *  differs from oldcomp;  use newcomp data: */
		if( ndiff != newcomp->natom ) {
			ndiff = newcomp->natom;
			if( diff == NULL ) {
				diff = (int *)malloc( ndiff*sizeof(int) );
			} else {
				diff = (int *)realloc( diff, ndiff*sizeof(int));
			}
			if( diff == NULL ) {
				fprintf( stderr,
				 "!!! put_ct: could not malloc diff\n");
				return MMIO_ERR;
			}
		}

		/* natom will be the number of atom records written: */
		if( (natom=filldiff(diff,oldcomp,newcomp)) < 0 ) {
			fprintf( stderr, "!!! put_ct: filldiff() fails\n" );
			return MMIO_ERR;
		}
		status = mmio_put_ct( iwrite, MMIO_COMPRESSED,
		 natom, newcomp->title);
		verbose( "put_ct: mmio_put_ct() returns %s\n",
		 mmio_return_code(status) );
		if( status == MMIO_ERR ) {
			fprintf( stderr, "!!! put_ct: mmio_put_ct() fails\n" );
			return status;
		}

		/* traverse the diff list, writing an atom record wherever
		 *  diff[] is true; note that since we are writing a
		 *  compressed CT, we can pass junk values to mmio_put_atom
		 *  for variables that are associated with connectivity: */
		for( idiff=0; idiffatom + idiff;
				status = mmio_put_atom(
				 iwrite, comp_atom->mmod_iatom,
				 0, 0,
				 NULL, NULL,
				 comp_atom->xyz, 0, 0,
				 '\0', comp_atom->color, 0,
				 '\0', NULL,
				 NULL );
				verbose(
				 "put_ct: mmio_put_atom() returns %s,"
				 " iatom= %d\n",
				 mmio_return_code(status), iatom );
				if( status == MMIO_ERR ) {
					fprintf( stderr,
					 "!!! put_ct: mmio_put_atom() fails\n");
					return status;
				}
			}
		}
	} else {
		/* write full CT, taking connectivity data from conn and
		 *  other data from oldcomp: */
		if( (natom=conn->natom) != oldcomp->natom ) {
			fprintf( stderr, "!!! put_ct: conn and comp have"
			 " different numbers of atoms\n" );
			return MMIO_ERR;
		}
		status = mmio_put_ct( iwrite, MMIO_FULL, natom, oldcomp->title);
		verbose( "put_ct: mmio_put_ct() returns %s\n",
		 mmio_return_code(status) );
		if( status == MMIO_ERR ) {
			fprintf( stderr, "!!! put_ct: mmio_put_ct() fails\n" );
			return status;
		}
		for( iatom=0; iatomatom + iatom;
			comp_atom = oldcomp->atom + iatom;
			status = mmio_put_atom( iwrite, comp_atom->mmod_iatom,
			 conn_atom->itype, conn_atom->nbond,
			 conn_atom->bond_atom, conn_atom->bond_order,
			 comp_atom->xyz, conn_atom->charge1, conn_atom->charge2,
			 conn_atom->chain, comp_atom->color, conn_atom->resnum,
			 conn_atom->resname1, conn_atom->resname4,
			 conn_atom->pdbname );
			verbose(
			 "put_ct: mmio_put_atom() returns %s, iatom= %d\n",
			 mmio_return_code(status), iatom );
			if( status == MMIO_ERR ) {
				fprintf( stderr,
				 "!!! put_ct: mmio_put_atom() fails\n" );
				return status;
			}
		}
	}
	return MMIO_OK;
}
/*******************************************************************/
static int filldiff( int AR diff, struct ct_comptag *oldcomp,
 struct ct_comptag *newcomp )
 /* array diff[] has the same number of elements as newcomp has atoms.
  * put TRUE in diff in positions where old and new atoms differ in
  *  coords or color;  FALSE where they don't;
  * return the number of places in which there are differences;
  *  this is the number of atom lines in the compressed CT which
  *  the calling function will write.
  * oldcomp is expected to contain data on all the atoms in the CT,
  *  and thus must have at least as many atoms as newcomp: */
{
	int ndiff = 0;
	int iatom, natom = newcomp->natom;
	struct comp_atomtag *oldatom = oldcomp->atom;
	struct comp_atomtag *end_oldatom = oldcomp->atom + oldcomp->natom;
	struct comp_atomtag *newatom = newcomp->atom;

	/* traverse the newatom list: */
	for( iatom=0; iatommmod_iatom != newatom->mmod_iatom ) {
			oldatom += 1;
			if( oldatom >= end_oldatom ) {
				fprintf( stderr, "!!! filldiff: oldatom end"
				 " encountered\n" );
				return -1;
			}
		}
		/* now compare newatom and oldatom data: */
		if( memcmp(oldatom->xyz,newatom->xyz,3*sizeof(float))!=0
		 || oldatom->color != newatom->color ) {
			diff[ iatom ] = TRUE;
			ndiff += 1;
		} else {
			diff[ iatom ] = FALSE;
		}
	}
	return ndiff;
}
/*******************************************************************/
Modified: Tue May 28 19:11:58 1996 GMT
Page accessed 4748 times since Sat Apr 17 21:57:58 1999 GMT