/aquad/src/AquaFuncts.c
C | 2571 lines | 1698 code | 345 blank | 528 comment | 396 complexity | 1a393f9bbaef0351b360e5314d264c18 MD5 | raw file
- /*
- ** File: AquaFuncts.c
- **
- ** This file is part of the AQUA program suite, which is being
- ** developed in the NMR Spectroscopy Department, Bijvoet Center for
- ** Biomolecular Research, Utrecht.
- **
- **
- ** AUTHOR(S) : Ton Rullmann
- ** VERSION : 1.5.0
- ** DATE : Dec 18, 1997
- **
- ** This file contains miscellaneous functions dealing with
- ** among other things Aqua data types.
- */
- #include <stdio.h>
- #include <string.h>
- #include <stdlib.h>
- #include <time.h>
- #include "Compiler.h"
- #include "General.h"
- #include "AquaMacros.h"
- #include "AquaTypes.h"
- #include "AquaData.h"
- /**************** function prototypes ****************/
- #include "AquaFuncts_cv.h"
- #include "AquaFuncts_io.h"
- #include "AquaStrucset.h"
- #include "AquaFiles.h"
- #include "MenuFuncts.h"
- #include "Functs.h"
- #include "AquaFuncts.h"
- static Atom *GuessAtom( Atom *at );
- static int CmpResidOfAtomptr( const void *a, const void *b );
- static Mol *ExpandAtoms( Atom *atoms, const int at_count, Mol *molptr );
- static Mol *MakeDefsMol( Atom *new_atoms, const int new_count, const Boolean makeres );
- /**************** globals ****************/
- extern Boolean iaf; /* import from Functs */
- extern int terseness; /* import from Functs */
- extern Card msg_string; /* import from Functs */
- char *project_name; /* export to main program and libs */
- #define max_MATCH 10
- static Atomptr MatchList[max_MATCH];
- static int find_atom_err; /* flag FindAtom2 -> WriteNotFoundAtom */
- /* BEGIN_EXTH */
- /**************** Check_General_Options ****************/
- /*
- ** Check the validity of the values of some general options against
- ** the allowed range
- */
- Boolean Check_General_Options ( void )
- {
-
- Boolean check = true;
- /* Terseness value [0-9] */
- if ( ( GetOptionIVal( opt_terseness ) < 0 ) ||
- ( GetOptionIVal( opt_terseness ) > 9 ) ) {
- check = false;
- printf("Option opt_terseness (%d) is outside it's bounds\n",GetOptionIVal(opt_terseness));
- }
-
- return ( check );
- }
- /**************** ProjectInit ****************/
- /*
- ** Get the project name from the environment variable AQUAPROJECT.
- */
- void ProjectInit( void )
- {
- project_name = GetVar( "AQUAPROJECT", false );
- if ( !project_name )
- {
- puts( " * Define a project name first!" );
- exit( 1 );
- }
- if ( terseness < 6 ) {
- printf( "Project name: %s\n", project_name );
- }
- return;
- }
- /**************** AquaDataFile ****************/
- /*
- ** Open Aqua data file <file> in the AQUADATADIR directory.
- */
- FILE *AquaDataFile( char *file )
- {
- FILE *fp;
- Filnam filnam;
- char *dir;
- /* the total length of the string should be checked against size_FILNAM */
- dir = GetVar( "AQUADATADIR", true );
- strcpy( filnam, dir );
- strcat( filnam, DIRSEP );
- strcat( filnam, file );
- if ( terseness < 6 ) {
- printf( "Now opening: %s\n", filnam );
- }
-
- if ( NOTOPENr( fp, filnam ) )
- {
- printf( " * Aqua data-file <%s> could not be opened\n", filnam );
- if ( !iaf )
- exit( 1 );
- return( NULL );
- }
- return( fp );
- }
- /**************** AquaFileName ****************/
- /*
- ** Generate an Aqua file name.
- ** The string is stored in <filnam> (which should contain at least
- ** size_FILNAM characters), and a pointer to it is returned.
- */
- char *AquaFileName( const Qtyp type, Typnam ident, char *filnam )
- {
- /* the total length of the string should be checked against size_FILNAM */
- strcpy( filnam, project_name );
- strcat( filnam, "_" );
- strcat( filnam, QtypSel( type ) );
- strcat( filnam, "_" );
- strcat( filnam, ident );
- strcat( filnam, FiltypExt( QtypFil( type ) ) );
-
- return( filnam );
- }
- /**************** OpenAquaFile ****************/
- /*
- ** Generate the name of an Aqua file of Qtyp <type> by inputting the "ident"
- ** (if <reuse> is false) or using the last input ident (if <reuse> is true)
- ** and calling AquaFileName,
- ** and open the file in mode <mode>.
- ** The return value is the file pointer, or null.
- ** <filnam> is the address of a string that receives the file name; it should
- ** contain at least size_FILNAM characters.
- */
- FILE *OpenAquaFile( const Qtyp type, char *mode, char *filnam, const Boolean reuse )
- {
- FILE *fp;
- Menu *oldmenu = NULL;
- static Typnam ident;
- char *filnam_temp;
- if ( !reuse )
- {
- oldmenu = SetOneMenu( "Give file identification string:" );
- MReadWord( ident );
- }
-
- filnam_temp = AquaFileName( type, ident, filnam );
-
- if ( terseness < 6 ) {
- if ( *mode == 'a' ) {
- printf( "Now (re)opening: %s\n", filnam_temp);
- } else {
- printf( "Now opening: %s\n", filnam_temp );
- }
- }
-
- if ( ( fp = fopen( filnam, mode ) ) == NULL )
- {
- printf( " * File <%s> could not be opened in mode <%s>\n", filnam, mode );
- if ( !iaf )
- exit( 1 );
- return( NULL );
- }
- if ( !reuse )
- ResetMenu( oldmenu );
- return( fp );
- }
- /**************** ReadAquaIDLines ****************/
- /*
- ** Read the standard Aqua identification lines.
- ** Returns true if successful, otherwise false.
- */
- Boolean ReadAquaIDLines( FILE *fp )
- {
- /***
- Card version;
- if ( fscanf( fp, "$ Aqua Version %s", version ) < 1 )
- {
- puts( " * This is not an Aqua file" );
- return( false );
- }
- if ( NOTEQ( version, aqua_version ) )
- {
- puts( " * This file was created by an old Aqua version" );
- return( false );
- }
- ***/
- float version;
- if ( fscanf( fp, "$ Aqua Version %f", &version ) < 1 )
- {
- puts( " * This is not an Aqua file" );
- return( false );
- }
- if ( version < aqua_versnum_comp )
- {
- puts( " * This file was created by an old Aqua version" );
- return( false );
- }
- else if ( version > aqua_versnum )
- {
- puts( " * This file was created by a future Aqua version" );
- exit( 1 );
- }
- FlushLine( fp );
- FlushLine( fp );
- return( true );
- }
- /**************** WriteAquaIDLines ****************/
- /*
- ** Write the standard Aqua identification lines
- */
- void WriteAquaIDLines( FILE *fp, const Qtyp qtype, Filnam filnam )
- {
- time_t tim;
- tim = time( NULL );
- fprintf( fp, "$ Aqua Version %s DATA: %s\n", aqua_version,
- QtypTxt( qtype ) );
- fprintf( fp, "$ FILE: %s CREATED: %s", filnam, ctime( &tim ) );
- }
- /**************** ProcessAquaLine ****************/
- /*
- ** Read until a non-blank line and process as standard Aqua header line.
- ** Format assumed:
- ** $ DESCRIPTOR: DATA
- ** The line is input stored at <line>.
- ** DESCRIPTOR and DATA are stored at <descr> and <inp>.
- ** These 3 first arguments should be pointers to non-overlapping
- ** character arrays of at least <size> bytes.
- ** The function returns NULL if an error occurs, otherwise <line>.
- */
- char *ProcessAquaLine( char *line, char *descr, char *inp,
- const int size, FILE *fp )
- {
- char *colptr;
- do
- if ( !ReadLine( line, size, fp ) )
- {
- puts( " * Error while reading" );
- return( NULL );
- }
- while
- ( IsEmptyStr( line ) );
- colptr = strchr( line, ':' );
- if ( *line != '$' || !colptr )
- {
- puts( " * Unexpected line encountered:" );
- puts( line );
- return( NULL );
- }
- strcpy( descr, FirstNonBlank( line+1 ) );
- strcpy( inp, FirstNonBlank( colptr+1 ) );
- if ( OptionTrue( opt_prtaq ) )
- printf( " read: %s\n", line );
- return( line );
- }
- /**************** QtypSel ****************/
- /*
- ** Return file selector string belonging to aqua type <type>
- */
- char *QtypSel( const Qtyp type )
- {
- unsigned int i;
-
- for ( i=0; i<size_aquarec; i++)
- if ( type == aquarec[i].type )
- return( aquarec[i].file_selector );
-
- return( aquarec[0].file_selector );
- }
- /**************** QtypTxt ****************/
- /*
- ** Return type text belonging to aqua type <type>
- */
- char *QtypTxt( const Qtyp type )
- {
- unsigned int i;
-
- for ( i=0; i<size_aquarec; i++)
- if ( type == aquarec[i].type )
- return( aquarec[i].type_text );
-
- return( aquarec[0].type_text );
- }
- /**************** QtypFil ****************/
- /*
- ** Return file type belonging to aqua type <type>
- */
- Filtyp QtypFil( const Qtyp type )
- {
- unsigned int i;
- for ( i=0; i<size_aquarec; i++)
- if ( type == aquarec[i].type )
- return( aquarec[i].file_type );
-
- return( aquarec[0].file_type );
- }
- /**************** FiltypNam ****************/
- /*
- ** Return name belonging to file type <type>
- */
- char *FiltypNam( const Filtyp type )
- {
- unsigned int i;
-
- for ( i=0; i<size_filerec; i++)
- if ( type == filerec[i].type )
- return( filerec[i].type_name );
-
- return( filerec[0].type_name );
- }
- /**************** FiltypLIBNam ****************/
- /*
- ** Return library name belonging to file type <type>
- */
- char *FiltypLIBNam( const Filtyp type )
- {
- switch ( type )
- {
- case pdb_type:
- return( "pdb" );
- case pdbn_type:
- return( "pdbn" );
- case pdbx_type:
- return( "pdbx" );
- case pdbm_type:
- return( "pdbm" );
- case biosym_type:
- return( "biosym" );
- default:
- return( "unknown" );
- }
- }
- /**************** FiltypExt ****************/
- /*
- ** Return file name extension belonging to file type <type>
- */
- char *FiltypExt( const Filtyp type )
- {
- unsigned int i;
-
- for ( i=0; i<size_filerec; i++)
- if ( type == filerec[i].type )
- return( filerec[i].extension );
-
- return( filerec[0].extension );
- }
- /**************** ContactRangeCount ****************/
- /*
- ** Return number of ContactRange types defined.
- */
- int ContactRangeCount( void )
- {
- return( size_contactrec );
- }
- /**************** ContactRangeNam ****************/
- /*
- ** Return name belonging to ContactRange type <type>.
- */
- char *ContactRangeNam( const ContactRange type )
- {
- unsigned int i;
-
- for ( i=0; i<size_contactrec; i++)
- if ( type == contactrec[i].type )
- return( contactrec[i].type_name );
-
- return( contactrec[0].type_name );
- }
- /**************** ContactRangeStarNam ****************/
- /*
- ** Return NMR-STAR name belonging to ContactRange type <type>.
- */
- char *ContactRangeStarNam( const ContactRange type )
- {
- unsigned int i;
-
- for ( i=0; i<size_contactrec; i++)
- if ( type == contactrec[i].type )
- return( contactrec[i].star_type_name );
-
- return( contactrec[0].star_type_name );
- }
- /**************** ContactRangeIndex ****************/
- /*
- ** Return index belonging to ContactRange type <type>.
- */
- int ContactRangeIndex( const ContactRange type )
- {
- unsigned int i;
-
- for ( i=0; i<size_contactrec; i++)
- if ( type == contactrec[i].type )
- return( i );
-
- return( -1 );
- }
- /**************** MapAtomsInSet ****************/
- /*
- ** For a strucset map the 'external' atom names to the
- ** 'internal' (Aqua) ones, or vice versa (determined by <how>).
- ** If needed, struc_mol's are created and a loop over them is performed.
- ** Cf. MapAtomsInDRestraints.
- */
- void MapAtomsInSet( Strucset *set, const int how, char *lib_nam )
- {
- char *fct_name = "MapAtomsInSet";
- Struc *strucptr;
- Mol *molptr = set->molecule;
- Filnam predef_nam;
- sprintf( msg_string, "preparing to map atom names in strucset <%s>",
- set->name );
- if ( terseness < 6 ) {
- STATUSm;
- }
- if ( set->mapped )
- {
- if ( terseness < 6 ) {
- STATUS( "atom names already mapped" );
- }
- return;
- }
- /* predefine the library name for use with the "whatif" option */
- strcpy( predef_nam, project_name );
- strcat( predef_nam, FiltypExt( atlib_type ) );
- /* DO THIS HERE ? >>>>>>>>>>>>>>>>>>>>>>>> */
- if ( EQn( lib_nam, "whatif", 6 ) && set->struc_count > 1 )
- set->onemap = false;
- /* use the unique molecule */
- if ( set->onemap )
- MapAtomsInMol( molptr, how, lib_nam, predef_nam );
- /* or use the model-dependent molecules */
- /* THE MAP IN set->molptr LOSES ITS MEANING */
- else
- {
- strucptr = set->fst_struc;
- while ( strucptr )
- {
- if ( !strucptr->struc_mol ) /* create if necessary */
- {
- strucptr->struc_mol = MolCreate( set->atom_count );
- MolCopy( strucptr->struc_mol, molptr );
- AddAtomSToMol( strucptr->struc_mol, molptr->atoms, molptr->atom_count,
- set->atom_count );
- }
- if ( set->struc_count != 1 ) /* name as given by whnamchk */
- sprintf( predef_nam, "%s_%03i%s", project_name,
- strucptr->struc_num, FiltypExt( atlib_type ) );
- MapAtomsInMol( strucptr->struc_mol, how, lib_nam, predef_nam );
- strucptr = strucptr->next_struc;
- }
- }
- set->mapped = true;
- }
- /**************** MapAtomsInMol ****************/
- /*
- ** Map the 'external' atom names to the 'internal' (Aqua) ones,
- ** or vice versa (determined by <how>).
- ** Cf. MapAtomsIn1DR.
- ** The library name <lib_nam> is derived from the coordinate file type,
- ** if it is "". Names are just copied if <lib_nam> is "copy" or "aqua".
- **
- ** Name guessing - as implemented in MapAtoms->FindAtom - is not performed.
- */
- void MapAtomsInMol( Mol *molptr, const int how, char *lib_nam_orig, Filnam predef_nam )
- {
- char *fct_name = "MapAtomsInMol";
- Mol *atmap_mol = NULL, *atmap_mol2 = NULL;
- FctFindAt Find = NULL;
- Strucset tmpSet;
- Card lib_nam;
- char *lib_nam1 = lib_nam, *lib_nam2;
- /* make a list of atom name pairs from a library file */
- /* the library name is derived from the coordinate file type, if it is "" */
- /* if the file type is "aqua", names will just be copied */
- /* was a combination of libraries specified? if so, split the name on the + character */
- strcpy( lib_nam, lib_nam_orig );
- if ( ( lib_nam2 = strchr( lib_nam, '+' ) ))
- {
- *lib_nam2 = '\0'; /* this sets lib_nam1 to a correct string */
- lib_nam2++;
- }
-
- /*
- printf( "first lib: <%s>\n", lib_nam1 );
- if ( lib_nam2 )
- printf( "second lib: <%s>\n", lib_nam2 );
- else
- puts( "no second lib" );
- */
- if ( NOTEQ( lib_nam1, "copy" ) && NOTEQ( lib_nam1, "aqua" ) )
- {
- if ( EQ( lib_nam1, "" ) )
- lib_nam1 = FiltypLIBNam( molptr->inp_filtyp );
- atmap_mol = ProcessAtomLIB( molptr, lib_nam1, predef_nam );
- }
- /* second library - handling same as in MapAtomsIn1DR */
- if ( lib_nam2 && atmap_mol )
- {
- if ( IsEmptyStr( lib_nam2 ) )
- lib_nam2 = FiltypLIBNam( molptr->inp_filtyp );
- /* lib_nam2 = GetLIBName( item_NONE ); */
- if ( ( atmap_mol2 = ProcessAtomLIB( molptr, lib_nam2, predef_nam ) ) )
- {
- StrucsetInit( &tmpSet, "allDefs", false );
- StrucsetCreate( &tmpSet, atmap_mol->atom_count + atmap_mol2->atom_count, 0 );
- AddAtomSToSet( &tmpSet, atmap_mol->atoms, atmap_mol->atom_count );
- AddAtomSToSet( &tmpSet, atmap_mol2->atoms, atmap_mol2->atom_count );
- DeleteMol( atmap_mol2 );
- if ( OptionTrue( opt_find2 ) )
- MakeResidues( &tmpSet, false ); /* cf. MakeDefsMol */
- atmap_mol = tmpSet.molecule;
- sprintf( msg_string, "%i atom name definitions active for present molecule",
- atmap_mol->atom_count );
- if ( terseness < 6 ) {
- STATUSm;
- }
- if ( OptionTrue( opt_prtlib ) )
- WriteFAtoms( stdout, atmap_mol->atoms, atmap_mol->atom_count );
- }
- }
- /* print what is going to happen */
- if ( !atmap_mol )
- {
- if ( terseness < 6 ) {
- if ( how == int_to_ext )
- STATUS( "copying internal to external atom names" );
- else
- STATUS( "copying external to internal atom names" );
- }
- }
- else
- {
- if ( how == int_to_ext )
- {
- if ( terseness < 6 ) {
- STATUS( "mapping internal to external atom names" );
- }
- if ( OptionTrue( opt_find2 ) )
- {
- if ( terseness < 6 ) {
- STATUS( "using routine FindAtom2" );
- }
- Find = FindAtom2;
- }
- else
- {
- if ( terseness < 6 ) {
- STATUS( "using routine FindAtom3" );
- }
- Find = FindAtom3;
- }
- }
- else
- {
- if ( terseness < 6 ) {
- STATUS( "mapping external to internal atom names" );
- }
- if ( OptionTrue( opt_find2 ) )
- {
- if ( terseness < 6 ) {
- STATUS( "using routine FindAtom2" );
- }
- Find = FindAtom2;
- }
- else
- {
- if ( terseness < 6 ) {
- STATUS( "using routine FindAtom3" );
- }
- Find = FindAtom3;
- }
- }
- }
- /* now do the mapping or copying */
- MapAtoms( molptr->atoms, molptr->atom_count, atmap_mol, how, false, Find );
- /* free the list again */
- if ( atmap_mol )
- DeleteMol( atmap_mol );
- }
- /**************** MapAtoms ****************/
- /*
- ** In <atoms> map the 'external' atom names to the 'internal' (Aqua)
- ** ones, or v.v., according to the names of the atoms appearing in <atmap_mol>.
- */
- void MapAtoms( Atom *atoms, int atom_count, Mol *atmap_mol,
- const int how, const Boolean guess, FctFindAt FindAtom )
- {
- Atom *atmap;
- Boolean err = false;
- int idx, find_which;
- if ( !atmap_mol )
- {
- /* if no list was made, just copy */
- if ( how == int_to_ext )
- {
- while ( atom_count-- )
- {
- memcpy( atoms->atom_ext, atoms->atom_int, size_ID );
- atoms++;
- }
- }
- else
- {
- while ( atom_count-- )
- {
- memcpy( atoms->atom_int, atoms->atom_ext, size_ID );
- atoms++;
- }
- }
- return;
- }
- /* if a list was made, use it to map */
- find_which = ( how == int_to_ext ? use_int : use_ext );
- while ( atom_count-- )
- {
- if ( !( atmap = FindAtom( &idx, atoms, atmap_mol, find_which, guess ) ) )
- {
- WriteNotFoundAtom( stdout, atoms, find_which, find_which );
- puts( "copying atom name instead" );
- if ( how == int_to_ext )
- memcpy( atoms->atom_ext, atoms->atom_int, size_ID );
- else
- memcpy( atoms->atom_int, atoms->atom_ext, size_ID );
- err = true;
- }
- else
- {
- if ( how == int_to_ext )
- memcpy( atoms->atom_ext, atmap->atom_ext, size_ID );
- else
- memcpy( atoms->atom_int, atmap->atom_int, size_ID );
- }
- atoms++;
- }
-
- if ( err )
- puts( " * One or more atoms could not be found in the AtomLIB" );
- }
- /**************** MapPtrAtoms ****************/
- /*
- ** In <atomptrs> map the 'external' atom names to the 'internal' (Aqua)
- ** ones, or v.v., according to the names of the atoms appearing in <atmap_mol>.
- */
- void MapPtrAtoms( Atomptr *atomptrs, int atom_count, Mol *atmap_mol,
- const int how, const Boolean guess, FctFindAt FindAtom )
- {
- Atom *atmap;
- Boolean err = false;
- int idx, find_which;
- if ( !atmap_mol )
- {
- /* if no list was made, just copy */
- if ( how == int_to_ext )
- {
- while ( atom_count-- )
- {
- memcpy( (*atomptrs)->atom_ext, (*atomptrs)->atom_int, size_ID );
- atomptrs++;
- }
- }
- else
- {
- while ( atom_count-- )
- {
- memcpy( (*atomptrs)->atom_int, (*atomptrs)->atom_ext, size_ID );
- atomptrs++;
- }
- }
- return;
- }
- /* if a list was made, use it to map */
- find_which = ( how == int_to_ext ? use_int : use_ext );
- while ( atom_count-- )
- {
- if ( !( atmap = FindAtom( &idx, *atomptrs, atmap_mol, find_which, guess ) ) )
- {
- WriteNotFoundAtom( stdout, *atomptrs, find_which, find_which );
- puts( "copying atom name instead" );
- if ( how == int_to_ext )
- memcpy( (*atomptrs)->atom_ext, (*atomptrs)->atom_int, size_ID );
- else
- memcpy( (*atomptrs)->atom_int, (*atomptrs)->atom_ext, size_ID );
- err = true;
- }
- else
- {
- if ( how == int_to_ext )
- memcpy( (*atomptrs)->atom_ext, atmap->atom_ext, size_ID );
- else
- memcpy( (*atomptrs)->atom_int, atmap->atom_int, size_ID );
- }
- atomptrs++;
- }
-
- if ( err )
- puts( " * One or more atoms could not be found in the AtomLIB" );
- }
- /**************** MakeResidues ****************/
- /*
- **
- */
- void MakeResidues( Strucset *set, Boolean sort )
- {
- char *fct_name = "MakeResidues";
- Mol *molptr;
- Atom *atoms;
- Atomptr *atomptrs, *list;
- Residue *residues;
- int *atomidxs, *list2;
- int i, j, idx, res_count, atom_count;
- char *res_id;
- molptr = set->molecule;
- atom_count = molptr->atom_count;
- /* make a list of pointers to the atoms */
- list = CREATE_NEW( Atomptr, atom_count );
- if ( !list )
- NOMEM( "atomptrs" );
- atomptrs = list;
- atoms = molptr->atoms;
- for ( i=0; i < atom_count; i++ )
- *atomptrs++ = atoms++;
- /* allocate a list of atom indices */
- list2 = CREATE_NEW( int, atom_count );
- if ( !list2 )
- NOMEM( "atomidxs" );
- /* make a unique residue ID name */
- atoms = molptr->atoms;
- for ( i=0; i < atom_count; i++ )
- MakeUniqResid( atoms++ );
- /* order the list on unique residue IDs */
- if ( sort )
- qsort( list, (size_t)atom_count, sizeof( Atomptr ), CmpResidOfAtomptr );
- /* now count the number of residues and set the new, sequential residue number */
- /* fill the atom index array */
- atomidxs = list2;
- atomptrs = list;
- atoms = *atomptrs;
- res_id = atoms->res_id;
- res_count = 1;
- for ( i=0; i < atom_count; i++ )
- {
- if ( NOTEQ( atoms->res_id, res_id ) )
- {
- res_id = atoms->res_id;
- res_count++;
- }
- atoms->res_num_new = res_count; /* sequential number, starting at 1 */
- for ( j=0; j < atom_count; j++ )
- {
- if ( atoms == &molptr->atoms[j] )
- {
- *atomidxs++ = j;
- break;
- }
- }
- atoms = *(++atomptrs);
- }
- /* make the residue list */
- sprintf( msg_string, "allocating %i residues in molecule of strucset <%s>",
- res_count, set->name );
- if ( terseness < 6 ) {
- STATUSm;
- }
- residues = CREATE_NEW( Residue, res_count );
-
- molptr->residues = residues;
- molptr->res_count = res_count;
- atomidxs = list2;
- atomptrs = list;
- atoms = *atomptrs;
- res_id = atoms->res_id;
- residues->res_num_new = atoms->res_num_new;
- residues->res_num_orig = atoms->res_num_orig;
- memcpy( residues->chain, atoms->chain, size_ID );
- memcpy( residues->res_nam, atoms->res_nam, size_ID );
- strcpy( residues->res_id, atoms->res_id );
- residues->res_insert = atoms->res_insert;
- residues->atomptrs = atomptrs;
- residues->atomidxs = atomidxs;
- for ( i=0; i < atom_count; i++ )
- {
- if ( NOTEQ( atoms->res_id, res_id ) )
- {
- residues++;
- res_id = atoms->res_id;
- residues->res_num_new = atoms->res_num_new;
- residues->res_num_orig = atoms->res_num_orig;
- memcpy( residues->chain, atoms->chain, size_ID );
- memcpy( residues->res_nam, atoms->res_nam, size_ID );
- strcpy( residues->res_id, atoms->res_id );
- residues->res_insert = atoms->res_insert;
- residues->atomptrs = atomptrs;
- residues->atomidxs = atomidxs;
- }
- residues->atom_count++;
- atoms = *(++atomptrs);
- ++atomidxs;
- }
- /* check */
- atoms = molptr->atoms;
- atomidxs = list2;
- atomptrs = list;
- for ( i=0; i < atom_count; i++ )
- {
- idx = *atomidxs++;
- if ( ( atoms + idx ) != *atomptrs++ )
- {
- ERROR( " * Atom index incorrect - code error" );
- exit( 1 );
- }
- }
- /* print */
- if ( OptionTrue( opt_prtres ) )
- WriteResidues( stdout, molptr->residues, molptr->res_count );
- }
- /**************** MakeUniqResid ****************/
- /*
- **
- */
- void MakeUniqResid( Atom *atom )
- {
- memcpy( atom->res_id, atom->chain, size_ID );
- sprintf( atom->res_id+size_ID, "%05d%c",
- atom->res_num_orig, atom->res_insert );
- memcpy( atom->res_id+size_ID+6, atom->res_nam, size_ID );
- *(atom->res_id+2*size_ID+6) = '\0';
- /*
- printf( "%s ", atom->res_id );
- PRINTID( atom->res_nam );
- putchar( '\n' );
- */
- }
- /**************** ReadUniqResid ****************/
- /*
- **
- */
- void ReadUniqResid( Atom *atom )
- {
- memcpy( atom->chain, atom->res_id, size_ID );
- memcpy( atom->res_nam, atom->res_id+size_ID+6, size_ID );
- if ( sscanf( atom->res_id+size_ID, "%5d%c",
- &atom->res_num_orig, &atom->res_insert ) < 2 )
- {
- printf( " * Cannot read from resid <%s>\n", atom->res_id );
- printf( " * Atom: " );
- WriteAtom( stdout, atom, use_both );
- putchar( '\n' );
- exit( 1 );
- }
- }
- /**************** WriteResidues ****************/
- /*
- **
- */
- void WriteResidues( FILE *fp, Residue *residues, const int count )
- {
- int i = 0;
- while ( i++ < count )
- {
- WriteResidue( fp, residues );
- fprintf( fp, "# atoms: %2i", residues->atom_count );
- fputs( ". first at: ", fp );
- WriteAtom( fp, *(residues->atomptrs), use_ext );
- fputs( ". last at: ", fp );
- WriteAtom( fp, *((residues->atomptrs)+residues->atom_count-1), use_ext );
- #ifdef DEBUG
- for ( j=0; j<residues->atom_count; j++ )
- fprintf( fp, " %i", residues->atomidxs[j] );
- #endif
- fputc( '\n', fp);
- residues++;
- }
- }
- /**************** CoordZero ****************/
- /*
- **
- */
- void CoordZero( Coord *v )
- {
- v->x = 0.F;
- v->y = 0.F;
- v->z = 0.F;
- }
- /**************** CoordAdd ****************/
- /*
- **
- */
- void CoordAdd( Coord *v, const Coord c )
- {
- v->x += c.x;
- v->y += c.y;
- v->z += c.z;
- }
- /**************** CoordDiv ****************/
- /*
- **
- */
- void CoordDiv( Coord *v, const int n )
- {
- v->x /= n;
- v->y /= n;
- v->z /= n;
- }
- /**************** FillAtom ****************/
- /*
- ** Fill an atom with data. Only one of the two names is set.
- **
- ** NOTE: in principle the new, sequential residue number could be found
- ** by going through the residue list, but for now comparisons will
- ** be based on the residue ID (cf. FindAtom).
- */
- Atom FillAtom( Resid res_id, ID res_nam, ID atom_nam, const int which )
- {
- /* char *fct_name = "FillAtom"; */
- Atom at;
- InitAtom( &at );
- if ( which == use_int )
- memcpy( at.atom_int, atom_nam, size_ID );
- if ( which == use_ext )
- memcpy( at.atom_ext, atom_nam, size_ID );
- memcpy( at.res_nam, res_nam, size_ID );
- strcpy( at.res_id, res_id );
- ReadUniqResid( &at );
- /*--at.res_num = at.res_num_orig; */
- /*
- STATUS("filling:");
- printf("%s ", res_id);
- PRINTIDv(res_nam);
- printf(" %i\n", at.res_num_orig);
- */
- return( at );
- }
- /**************** FindAtom1 ****************/
- /*
- ** Try to find an atom <at> in the external or internal (depending on <which>)
- ** atom name list of a molecule <molptr>.
- ** Residues are matched on 'res_id'.
- ** If <guess>, some name mappings are guessed (only for external names).
- ** The list index is returned in <idx>, and the corresponding atom in
- ** the function value.
- */
- Atom *FindAtom1( int *idx, Atom *at, Mol *molptr, const int which, const Boolean guess )
- {
- /* char *fct_name = "FindAtom1"; */
- int i = 0, count;
- /* int nmatch; */
- Atom *atoms, *match, *try_variable;
-
- atoms = molptr->atoms;
- count = molptr->atom_count;
- /* loop over all atoms and compare */
- #ifdef DEBUG
- STATUS( "trying to find:" );
- printf( " " );
- WriteAtom( stdout, at, use_both );
- putchar( '\n' );
- #endif
- while ( i < count )
- {
- #ifdef DEBUG
- printf( "looking at: " );
- WriteAtom( stdout, atoms, use_both );
- #endif
- /* match residue ID's */
- if ( EQ( atoms->res_id, at->res_id ) )
- {
- /* match atom names */
- if ( which == use_ext )
- {
- if ( EQID( atoms->atom_ext, at->atom_ext ) )
- break;
- }
- else if ( which == use_int )
- {
- if ( EQID( atoms->atom_int, at->atom_int ) )
- break;
- }
- }
- atoms++;
- i++;
- }
- *idx = i;
- /* found */
- if ( i < count )
- return( atoms );
- /* not found */
- if ( which == use_ext && guess )
- {
- /* derive a closely related name */
- if ( ( try_variable = GuessAtom( at ) ) == NULL )
- {
- /* find wildcard matches */
- #ifdef DEBUG
- FindAtomsWC( &nmatch, at, molptr, use_ext );
- #endif
- return( NULL );
- }
-
- /* try again */
- #ifdef DEBUG
- STATUS( "trying again with:" );
- printf( " " );
- WriteAtom( stdout, try_variable, use_both);
- putchar( '\n' );
- #endif
- match = FindAtom1( idx, try_variable, molptr, use_ext, false );
- if ( match )
- {
- printf( " External atom name remapped for: " );
- WriteAtom( stdout, at, use_ext );
- printf( " on: " );
- PRINTIDv( try_variable->atom_ext );
- putchar( '\n' );
- return( match );
- }
- else
- {
- /* find wildcard matches */
- #ifdef DEBUG
- FindAtomsWC( &nmatch, at, molptr, use_ext );
- #endif
- return( NULL );
- }
- }
- return( NULL );
- }
- /**************** FindAtom2 ****************/
- /*
- ** New version of FindAtom1.
- ** - Faster: makes use of residue list.
- ** - Correct residue names where needed.
- **
- ** Try to find an atom <at> in the external or internal (depending on <which>)
- ** atom name list of a molecule <molptr>.
- ** Residues are matched on 'res_id'.
- ** If <guess>, some name mappings are guessed (only for external names).
- ** The list index is returned in <idx>, and the corresponding atom in
- ** the function value.
- **
- ** If the residue list is scrambled (i.e. atoms belonging to one residue are
- ** scattered across more than one entry in molptr->residues) the routine
- ** fails (cf. comment in ExpandAtoms->MakeDefsMol).
- */
- Atom *FindAtom2( int *idx, Atom *theAtom, Mol *molptr, const int which,
- const Boolean guess )
- {
- /* char *fct_name = "FindAtom2"; */
- Atom *at = theAtom, *match, mod, *try_variable;
- Atomptr *atomptrs;
- Residue *res;
- char *res_nam, *lastcharptr;
- int *atomidxs, res_idx, count;
- /* int nmatch; */
-
- #ifdef DEBUG
- STATUS( "trying to find:" );
- printf( " " );
- WriteAtom( stdout, at, use_both );
- putchar( '\n' );
- #endif
- find_atom_err = 0;
- /* find the residue; map Arg+ onto Arg etc. */
- /*
- (The mapping is necessary because the expanded AtomLIB definitions will
- not contain Arg+ etc. (even though they exist in the unexpanded AtomLIB)
- if the set (i.e. the PDB file) does not contain Arg+ etc.; the mapping
- is only necessary to find the residue in the expanded AtomLIB.)
- */
- res_nam = at->res_nam;
- res_idx = FindResidueIndexFromId( at->res_id, molptr );
- if ( res_idx == -1 )
- {
- lastcharptr = LastNonBlank( at->res_nam, size_ID );
- if ( *lastcharptr == '+' || *lastcharptr == '-' )
- {
- printf( " Residue name remapped for: " );
- WriteAtom( stdout, at, which );
- memcpy( &mod, at, sizeof( Atom ) );
- at = &mod;
- *LastNonBlank( at->res_nam, size_ID ) = ' ';
- MakeUniqResid( at );
- res_idx = FindResidueIndexFromId( at->res_id, molptr );
- printf( " on: " );
- PRINTIDv( at->res_nam );
- putchar( '\n' );
- }
- }
- if ( res_idx == -1 )
- {
- find_atom_err = 1;
- return( NULL );
- }
- #ifdef DEBUG
- sprintf( msg_string, "found residue index %i", res_idx );
- STATUSm;
- #endif
- /* find the atom */
- res = &molptr->residues[res_idx];
- count = res->atom_count;
- atomptrs = res->atomptrs;
- atomidxs = res->atomidxs;
- while ( count )
- {
- #ifdef DEBUG
- printf( "looking at: " );
- WriteAtom( stdout, *atomptrs, use_both );
- printf( " index %i\n", *atomidxs );
- #endif
- if ( which == use_ext )
- {
- if ( EQID( (*atomptrs)->atom_ext, at->atom_ext ) )
- {
- *idx = *atomidxs;
- break;
- }
- }
- else if ( which == use_int )
- {
- if ( EQID( (*atomptrs)->atom_int, at->atom_int ) )
- {
- *idx = *atomidxs;
- break;
- }
- }
- atomptrs++;
- atomidxs++;
- count--;
- }
- /* found */
- if ( count )
- return( *atomptrs );
- /* not found */
- if ( which == use_ext && guess )
- {
- /* derive a closely related name */
- if ( ( try_variable = GuessAtom( at ) ) == NULL )
- {
- /* find wildcard matches */
- #ifdef DEBUG
- FindAtomsWC( &nmatch, at, molptr, use_ext );
- #endif
- return( NULL );
- }
-
- /* try again */
- #ifdef DEBUG
- STATUS( "trying again with:" );
- printf( " " );
- WriteAtom( stdout, try_variable, use_both);
- putchar( '\n' );
- #endif
- match = FindAtom2( idx, try_variable, molptr, use_ext, false );
- if ( match )
- {
- printf( " External atom name remapped for: " );
- WriteAtom( stdout, at, use_ext );
- printf( " on: " );
- PRINTIDv( try_variable->atom_ext );
- putchar( '\n' );
- return( match );
- }
- else
- {
- /* find wildcard matches */
- #ifdef DEBUG
- FindAtomsWC( &nmatch, at, molptr, use_ext );
- #endif
- return( NULL );
- }
- }
- return( NULL );
- }
- /**************** FindAtom3 ****************/
- /*
- ** New version of FindAtom1.
- ** Handles chain and residue number wildcards such as set by AtomLIB.
- ** This routine captures the matching of wildcards that was formerly
- ** handled by ExpandAtoms, and is meant to be used by MapAtoms/MapPtrAtoms
- ** only.
- ** Residue name variants (such as Arg+) must be explicitly present in the
- ** molecule stored at <molptr> (i.e. in the AtomLIB).
- ** In contrast to FindAtom2 the molecule does not have to contain a
- ** residue list.
- */
- Atom *FindAtom3( int *idx, Atom *at, Mol *molptr, const int which,
- const Boolean guess )
- {
- /* char *fct_name = "FindAtom3"; */
- int i = 0, count;
- /* int nmatch; */
- Atom *atoms, *match, *try_variable;
-
- atoms = molptr->atoms;
- count = molptr->atom_count;
- /* loop over all atoms and compare */
- #ifdef DEBUG
- STATUS( "trying to find:" );
- printf( " " );
- WriteAtom( stdout, at, use_both );
- putchar( '\n' );
- #endif
- while ( i < count )
- {
- #ifdef DEBUG
- printf( "looking at: " );
- WriteAtom( stdout, atoms, use_both );
- #endif
- /* match residue and chain */
- if ( NOTEQID( atoms->res_nam, at->res_nam ) )
- {
- atoms++;
- i++;
- continue;
- }
- if ( NOTEQID( atoms->chain, at->chain ) && atoms->chain[0] != '*' )
- {
- atoms++;
- i++;
- continue;
- }
- if ( atoms->res_num_orig == at->res_num_orig ||
- atoms->res_num_orig == wild_RESNUM )
- {
- /* match atom names */
- /* specific defs MUST come BEFORE wildcards in order to be found first */
- if ( which == use_ext )
- {
- if ( EQID( atoms->atom_ext, at->atom_ext ) )
- break;
- }
- else if ( which == use_int )
- {
- if ( EQID( atoms->atom_int, at->atom_int ) )
- break;
- }
- }
- atoms++;
- i++;
- }
- *idx = i;
- /* found */
- if ( i < count )
- return( atoms );
- /* not found */
- if ( which == use_ext && guess )
- {
- /* derive a closely related name */
- if ( ( try_variable = GuessAtom( at ) ) == NULL )
- {
- /* find wildcard matches */
- #ifdef DEBUG
- FindAtomsWC( &nmatch, at, molptr, use_ext );
- #endif
- return( NULL );
- }
-
- /* try again */
- #ifdef DEBUG
- STATUS( "trying again with:" );
- printf( " " );
- WriteAtom( stdout, try_variable, use_both);
- putchar( '\n' );
- #endif
- match = FindAtom3( idx, try_variable, molptr, use_ext, false );
- if ( match )
- {
- printf( " External atom name remapped for: " );
- WriteAtom( stdout, at, use_ext );
- printf( " on: " );
- PRINTIDv( try_variable->atom_ext );
- putchar( '\n' );
- return( match );
- }
- else
- {
- /* find wildcard matches */
- #ifdef DEBUG
- FindAtomsWC( &nmatch, at, molptr, use_ext );
- #endif
- return( NULL );
- }
- }
- return( NULL );
- }
- /**************** TestFindAtom2 ****************/
- /*
- **
- */
- void TestFindAtom2( Mol *molptr )
- {
- char *fct_name = "TestFindAtom2";
- int idx, count = molptr->atom_count;
- Atom *atoms = molptr->atoms, *at;
-
- STATUS( "entering" );
- while ( count-- )
- {
- if ( ! ( at = FindAtom2( &idx, atoms, molptr, use_ext, true ) ) )
- WriteNotFoundAtom( stdout, atoms, use_ext, use_both );
- printf( "%i ", idx );
- if ( at != molptr->atoms + idx )
- {
- puts( "error" );
- WriteAtom( stdout, at, use_both ); putchar( '\n' );
- WriteAtom( stdout, &atoms[idx], use_both ); putchar( '\n' );
- }
- else
- putchar( '\n' );
- atoms++;
- }
- STATUS( "leaving" );
- }
- /**************** FindAtomsWC ****************/
- /*
- ** Try to find an atom <at> in the external or internal (depending on <which>
- ** atom name list of a molecule <molptr>,
- ** matching UP TO a wildcard character # or *.
- ** Residues are matched on 'res_id'.
- ** The number of matches is returned in <nmatch>, and the address of an
- ** Atomptr array in the function value.
- */
- Atomptr *FindAtomsWC( int *nmatch, Atom *at, Mol *molptr, const int which )
- {
- /* char *fct_name = "FindAtomsWC"; */
- char *wc;
- size_t len = 0;
- char test[size_ID+1];
- Atom *atoms;
- Atomptr *match = MatchList;
- int i=0, n=0;
- int count;
-
- atoms = molptr->atoms;
- count = molptr->atom_count;
- *nmatch = 0;
- /* loop over all atoms and compare */
- if ( which == use_ext )
- {
- strncpy( test, at->atom_ext, size_ID );
- if ( ( wc = strpbrk( test, "#*" ) ) == NULL )
- return( NULL );
- len = (size_t)( wc - test );
- }
- else if ( which == use_int )
- {
- strncpy( test, at->atom_int, size_ID );
- if ( ( wc = strpbrk( test, "#*" ) ) == NULL )
- return( NULL );
- len = (size_t)( wc - test );
- }
- while ( i < count )
- {
- /* match residue ID's -- assume unordered list */
- if ( EQ( atoms->res_id, at->res_id ) )
- {
- /* match atom names */
- if ( which == use_ext )
- {
- if ( EQn( atoms->atom_ext, at->atom_ext, len ) )
- {
- if ( ++n > max_MATCH )
- {
- puts( " * ERROR: max_MATCH too small" );
- exit( 1 );
- }
- *match++ = atoms;
- }
- }
- else if ( which == use_int )
- {
- if ( EQn( atoms->atom_int, at->atom_int, len ) )
- {
- if ( ++n > max_MATCH )
- {
- puts( " * ERROR: max_MATCH too small" );
- exit( 1 );
- }
- *match++ = atoms;
- }
- }
- }
- atoms++;
- i++;
- }
- *nmatch = n;
- #ifdef DEBUG
- if ( n ) {
- STATUS( "found wildcard representations:" );
- for ( i=0; i<n; i++ )
- {
- printf( " " );
- WriteAtom( stdout, MatchList[i], use_both );
- putchar('\n');
- }
- }
- #endif
- if ( n )
- return( MatchList );
- else
- return( NULL );
- }
- /**************** WriteNotFoundAtom ****************/
- /*
- **
- */
- void WriteNotFoundAtom( FILE *fp, Atom *atom, const int find_which,
- const int print_which )
- {
- if ( find_atom_err == 1 )
- fputs( " Residue name not found: ", fp );
- else
- {
- if ( find_which == use_ext )
- fputs( " External atom name not found: ", fp );
- if ( find_which == use_int )
- fputs( " Internal atom name not found: ", fp );
- }
- WriteAtom( fp, atom, print_which );
- fputc( '\n', fp );
- }
- /**************** FindResidueIndexFromNum ****************/
- /*
- ** Returns the index in the residues-array of the residue
- ** for which the new residue number is <res_num>.
- */
- int FindResidueIndexFromNum( const int res_num, Mol *molptr )
- {
- int i = 0, count;
- Residue *residues;
-
- residues = molptr->residues;
- count = molptr->res_count;
- while ( i < count )
- {
- if ( res_num == residues->res_num_new )
- return( i );
- residues++;
- i++;
- }
- return( -1 );
- }
- /**************** FindResidueIndexFromOrigNum ****************/
- /*
- ** Returns the index in the residues-array of the residue
- ** for which the original residue number, chain id and insertion code are
- ** <res_num>, <chain> and <insert>.
- */
- int FindResidueIndexFromOrigNum( const int res_num, ID chain, char insert,
- Mol *molptr )
- {
- int i = 0, count;
- Residue *residues;
-
- residues = molptr->residues;
- count = molptr->res_count;
- while ( i < count )
- {
- if ( res_num == residues->res_num_orig &&
- EQID( chain, residues->chain ) && insert == residues->res_insert )
- return( i );
- residues++;
- i++;
- }
- return( -1 );
- }
- /**************** FindResidueIndexFromId ****************/
- /*
- ** Returns the residue index of the residue
- ** for which the ID is <res_id>, by searching the residue list.
- */
- int FindResidueIndexFromId( Resid res_id, Mol *molptr )
- {
- int i = 0, count;
- Residue *residues;
-
- residues = molptr->residues;
- count = molptr->res_count;
- while ( i < count )
- {
- if ( EQ( res_id, residues->res_id ) )
- return( i );
- residues++;
- i++;
- }
- return( -1 );
- }
- /**************** FindResidueNumFromId ****************/
- /*
- ** Returns the new residue number of the residue
- ** for which the ID is <res_id>, by searching the residue list.
- */
- int FindResidueNumFromId( Resid res_id, Mol *molptr )
- {
- int i = 0, count;
- Residue *residues;
-
- residues = molptr->residues;
- count = molptr->res_count;
- while ( i < count )
- {
- if ( EQ( res_id, residues->res_id ) )
- return( residues->res_num_new );
- residues++;
- i++;
- }
- return( -1 );
- }
- /**************** FindResidueNamFromId ****************/
- /*
- ** Returns the residue name of the residue for which the ID
- ** is <res_id> (not comparing the part which corresponds to the residue name
- ** of course...) by searching the residue list.
- */
- char *FindResidueNamFromId( Resid res_id, Mol *molptr )
- {
- int i = 0, count;
- Residue *residues;
- static char unk[size_ID+1] = "UNK ";
-
- residues = molptr->residues;
- count = molptr->res_count;
- while ( i < count )
- {
- if ( EQn( res_id, residues->res_id, size_ID+6 ) )
- return( residues->res_nam );
- residues++;
- i++;
- }
- return( unk );
- }
- /**************** InitAtomS ****************/
- /*
- **
- */
- void InitAtomS( Atom *atom, int count )
- {
- while ( count-- )
- InitAtom( atom++ );
- }
- /**************** InitAtom ****************/
- /*
- **
- */
- void InitAtom( Atom *atom )
- {
- atom->res_num_new = 0;
- atom->res_num_orig = 0;
- Blank( atom->chain, size_ID );
- Blank( atom->atom_int, size_ID );
- Blank( atom->atom_ext, size_ID );
- Blank( atom->atom_element_symbol, size_ID );
- Blank( atom->res_nam, size_ID );
- atom->res_insert = ' ';
- *(atom->res_id) = '\0';
- }
- /**************** InitPsatomS ****************/
- /*
- **
- */
- void InitPsatomS( Psatom *psatom, int count )
- {
- while ( count-- )
- InitPsatom( psatom++ );
- }
- /**************** InitPsatom ****************/
- /*
- **
- */
- void InitPsatom( Psatom *psatom )
- {
- psatom->pseudo_typ = 0;
- psatom->num_form = 0;
- psatom->res_num_orig = 0;
- Blank( psatom->chain, size_ID );
- Blank( psatom->res_nam, size_ID );
- Blank( psatom->pseudo_nam, size_ID );
- Blank( (char *)psatom->form_nam, size_ID*max_PSFORM );
- psatom->res_insert = ' ';
- *(psatom->res_id) = '\0';
- }
- /**************** GetLIBName ****************/
- /*
- ** Determine the library name for the AtomLIB to be used.
- ** If option <option> has been set (i.e. <option> is not item_NONE and
- ** the option value is not "none"), the routine uses the option value
- ** instead of presenting a menu.
- */
- char *GetLIBName( const int option )
- {
- Menu *oldmenu;
- int item;
- char *name, *optval;
- static MenuItem menulist[] =
- {
- { item_aqua, "aqua", none, on },
- { item_disman, "disman", none, on },
- { item_disgeo, "disgeo", none, on },
- { item_diana, "diana", none, on },
- { item_xplor, "xplor", none, on },
- { item_xplori, "xplor-inv", none, on },
- { item_biosym, "biosym", none, on },
- { item_whatif, "whatif", none, on },
- { item_local, "local", none, on },
- { item_user, "user", none, on },
- { item_copy, "copy", none, on },
- { item_END, "end", none, on } /* same as in AquaWhat, AquaPseudo and AquaHow */
- };
- static Menu menu =
- {
- "GetLIBName menu", "format", "Give atom name format:",
- sizeof(menulist)/sizeof(MenuItem)
- };
-
- oldmenu = SetMenu( &menu, menulist );
- if ( option != item_NONE )
- optval = GetOptionCVal( option ); /* use value of option */
- /* OK if option menu has not been swapped */
- else
- optval = "none";
- while ( true )
- {
- if ( EQ( optval, "none" ) )
- item = ReadMenuItem();
- else
- item = FindMenuItem( optval );
- switch ( item )
- {
- case item_aqua:
- case item_disman:
- case item_disgeo:
- case item_diana:
- case item_xplor:
- case item_xplori:
- case item_biosym:
- case item_whatif:
- case item_local:
- case item_user:
- case item_copy:
- name = GetMenuItemKey( item );
- ResetMenu( oldmenu );
- return( name );
- case item_END:
- case item_NULL:
- ResetMenu( oldmenu );
- return( NULL );
- default:
- ProcessMenuError( menu_fatal );
- break;
- }
- }
- }
- /**************** ProcessAtomLIB ****************/
- /*
- ** Read and process an AtomLIB i.e. a file with atom name mappings.
- */
- Mol *ProcessAtomLIB( Mol *molptr, char *lib_nam, Filnam predef )
- {
- char *fct_name = "ProcessAtomLIB";
- char *usrnam = NULL;
- Atom *atmap_list;
- Mol *atmap_mol;
- Filnam lib_filnam;
- FILE *lib_fp = NULL;
- Boolean done = false;
- Menu *oldmenu;
- int atmap_count;
- if ( EQ( lib_nam, "user" ) || EQ( lib_nam, "local" ) || EQ( lib_nam, "whatif" ) )
- /* user defined library in Aqua format */
- {
- oldmenu = SetOneMenu( "Give name of atom-library file (or - to end):" );
- if ( EQ( lib_nam, "local" ) )
- usrnam = "AtomLIB";
- else if ( EQ( lib_nam, "whatif" ) )
- usrnam = predef;
- else
- usrnam = GetOptionCVal( opt_atmlib ); /* use value of atomlib option */
- /* OK if option menu has not been swapped */
- while ( !done )
- {
- if ( !usrnam )
- {
- MReadWord( lib_filnam );
- if ( *lib_filnam == '-' || strlen( lib_filnam ) == 0 )
- {
- ResetMenu( oldmenu );
- return( NULL );
- }
- }
- else
- {
- strcpy( lib_filnam, usrnam );
- usrnam = NULL; /* to prevent infinite loop */
- }
- if ( terseness < 6 ) {
- printf( "Now opening: %s\n", lib_filnam );
- }
- if ( ! ( done = OPENr( lib_fp, lib_filnam ) ) )
- {
- printf( " * atom-library file <%s> could not be opened\n", lib_filnam );
- MError();
- }
- }
- ResetMenu( oldmenu );
- if ( terseness < 6 ) {
- STATUS( "reading user AtomLIB" );
- }
- }
- else
- /* standard Aqua library from Aqua data directory */
- {
- strcpy( lib_filnam, "AtomLIB-" );
- strcat( lib_filnam, lib_nam );
- lib_fp = AquaDataFile( lib_filnam );
- if ( !lib_fp )
- return( NULL );
- if ( terseness < 6 ) {
- STATUS( "reading standard AtomLIB" );
- }
- }
- atmap_list = ReadAtomLIB( lib_fp, &atmap_count );
- if ( !atmap_list || atmap_count == 0 )
- {
- puts( " * Atom name definition list empty" );
- return( NULL );
- }
- sprintf( msg_string, "%i atom name definitions read", atmap_count );
- if ( terseness < 6 ) {
- STATUSm;
- }
- if ( OptionTrue( opt_prtlib ) )
- WriteFAtoms( stdout, atmap_list, atmap_count );
- if ( !OptionTrue( opt_find2 ) )
- return( MakeDefsMol( atmap_list, atmap_count, false ) );
- /* expand the list from 'wildcard' residue numbers to individual residues */
- if ( terseness < 6 ) {
- STATUS( "generating definitions for present molecule" );
- }
- atmap_mol = ExpandAtoms( atmap_list, atmap_count, molptr );
- free( atmap_list );
- sprintf( msg_string, "%i atom name definitions active for present molecule",
- atmap_mol->atom_count );
- if ( terseness < 6 ) {
- STATUSm;
- }
- if ( OptionTrue( opt_prtlib ) )
- WriteFAtoms( stdout, atmap_mol->atoms, atmap_mol->atom_count );
- return( atmap_mol );
- }
- /**************** ProcessStereoLIB ****************/
- /*
- ** Read and process a StereoLIB i.e. a file describing the "equivalent
- ** atoms" in prochiral groups.
- ** This routine is an analogue of ProcessPseudo.
- ** The descriptions and data types are the same as for pseudo atoms.
- */
- Psatom *ProcessStereoLIB( int *stereo_count, Mol *molptr, char *lib )
- {
- char *usrnam = NULL, *fct_name = "ProcessStereoLIB";
- char *lib_nam;
- Psatom *stereo_list;
- Filnam lib_filnam;
- FILE *lib_fp = NULL;
- Boolean done = false;
- Menu *oldmenu;
- lib_nam = GetLIBName( opt_rstinp );
- if ( !lib_nam )
- return( NULL ); /* no processing done */
- if ( EQ( lib_nam, "user" ) )
- /* type "user": user defined library in Aqua format */
- {
- oldmenu = SetOneMenu( "Give name of stereo-library file:" );
- usrnam = GetOptionCVal( opt_sterlib ); /* use value of stereolib option */
- /* OK if option menu has not been swapped */
- while ( !done )
- {
- if ( !usrnam )
- MReadWord( lib_filnam );
- else
- {
- strcpy( lib_filnam, usrnam );
- usrnam = NULL; /* to prevent infinite loop */
- }
- if ( terseness < 6 ) {
- printf( "Now opening: %s\n", lib_filnam );
- }
- if ( ! ( done = OPENr( lib_fp, lib_filnam ) ) )
- {
- printf( " * stereo-library file <%s> could not be opened\n", lib_filnam );
- MError();
- }
- }
- ResetMenu( oldmenu );
- if ( terseness < 6 ) {
- STATUS( "Reading user StereoLIB" );
- }
- }
- else
- /* type "local", "diana", etc: standard Aqua library from Aqua data directory */
- {
- lib_fp = AquaDataFile( lib );
- if ( !lib_fp )
- exit( 1 );
- if ( terseness < 6 ) {
- STATUS( "Reading standard StereoLIB" );
- }
- }
- stereo_list = ReadPseudoLIB( lib_fp, stereo_count );
- if ( !stereo_list || *stereo_count == 0 )
- {
- puts( " * Stereo definition list empty" );
- return( NULL );
- }
- /*
- sprintf( msg_string, "%i stereo definitions read", *stereo_count );
- STATUSm;
- PrintPsatomS( stereo_list, *stereo_count );
- */
- /* expand the list from 'wildcard' residue numbers to individual residues */
- /* NECESSARY?? cf. notes about MapAtoms.... */
- if ( terseness < 6 ) {
- STATUS( "generating definitions for present molecule" );
- }
- stereo_list = ExpandPsatoms( stereo_list, stereo_count, molptr );
- sprintf( msg_string, "%i stereo definitions active for present molecule",
- *stereo_count );
- if ( terseness < 6 ) {
- STATUSm;
- }
-
- if ( OptionTrue( opt_prtlib ) && ( terseness < 9 ) ) {
- PrintPsatomS( stereo_list, *stereo_count );
- }
-
- return( stereo_list );
- }
- /**************** ExpandPsatoms ****************/
- /*
- ** Expand the Psatom list from 'wildcard' residue numbers to individual
- ** residues.
- **
- ** NOTE: res_num is treated as res_num_orig
- ** chain id not yet handled
- ** Handling is such that for ** wildcarded residues **
- ** occurrences in every chain are handled, and the residue info
- ** (including ID) is completed ** by copying from the residue list **.
- **
- ** Mechanism: see ExpandAtoms.
- */
- Psatom *ExpandPsatoms( Psatom *psatoms, int *pseudo_count, Mol *molptr )
- {
- char *fct_name = "ExpandPsatoms";
- int count, new_count = 0, res_count;
- Psatom *list, *new_list, *new_psatoms;
- Residue *resptr;
- /* determine how many atoms are in the molecule for each 'wildcarded' residue number */
- /* use same loop mechanism as in actual setup of list (cf. below) */
- res_count = molptr->res_count;
- resptr = molptr->residues;
- while ( res_count-- )
- {
- list = psatoms;
- count = *pseudo_count;
- while ( count-- )
- {
- if ( NOTEQID( resptr->res_nam, list->res_nam ) )
- {
- list++;
- continue;
- }
- if ( resptr->res_num_orig == list->res_num_orig ||
- list->res_num_orig == wild_RESNUM )
- {
- new_count++;
- }
- list++;
- }
- resptr++;
- }
- /* now setup the new list */
- new_psatoms = CREATE_NEW( Psatom, new_count );
- if ( !new_psatoms )
- NOMEM( "new_psatoms" );
- new_list = new_psatoms;
- res_count = molptr->res_count;
- resptr = molptr->residues;
- while ( res_count-- )
- {
- list = psatoms;
- count = *pseudo_count;
- while ( count-- )
- {
- if ( NOTEQID( resptr->res_nam, list->res_nam ) )
- {
- list++;
- continue;
- }
- if ( resptr->res_num_orig == list->res_num_orig ||
- list->res_num_orig == wild_RESNUM )
- /* specific defs MUST come BEFORE wildcards */
- {
- /* copy the atom description */
- memcpy( new_list, list, sizeof( Psatom ) );
- /* fill in the residue info */
- memcpy( new_list->chain, resptr->chain, size_ID );
- strcpy( new_list->res_id, resptr->res_id );
- new_list->res_insert = resptr->res_insert;
- new_list->res_num_orig = resptr->res_num_orig;
- new_list++;
- }
- list++;
- }
- resptr++;
- }
- /* destroy the old list, and return the new list and new count */
- free( psatoms );
- *pseudo_count = new_count;
- return( new_psatoms );
- }
- /**************** PrintPsatomS ****************/
- /*
- ** Print a list of Psatom's.
- */
- void PrintPsatomS( Psatom *psatoms, const int pseudo_count )
- {
- int count = 0;
- puts( "List of pseudo atom definitions:" );
- while ( count++ < pseudo_count )
- {
- printf( "%i: ", count );
- PrintPsatom( psatoms++ );
- }
- }
- /**************** PrintPsatom ****************/
- /*
- ** Print a Psatom.
- */
- void PrintPsatom( Psatom *psatoms )
- {
- int i;
- if ( LenChars( psatoms->chain, size_ID ) )
- {
- printf( "Chain " );
- PRINTIDv( psatoms->chain );
- putchar( ' ' );
- }
- PRINTIDv( psatoms->res_nam );
- if ( psatoms->res_num_orig == wild_RESNUM )
- printf( " * typ %i ", psatoms->pseudo_typ );
- else
- printf( " %i%c typ %i ", psatoms->res_num_orig, psatoms->res_insert,
- psatoms->pseudo_typ );
- PRINTIDv( psatoms->pseudo_nam );
- putchar( ' ' );
- i = 0;
- while ( i < psatoms->num_form )
- {
- PRINTIDv( psatoms->form_nam[i] );
- i++;
- putchar( ' ' );
- }
- putchar( '\n' );
- }
- /* END_EXTH ----------------------- only private functions follow */
- /**************** GuessAtom ****************/
- /*
- ** Returns modified version of <at> (change of external name)
- ** or NULL if no suitable mapping exists.
- */
- static Atom *GuessAtom( Atom *at )
- {
- static Atom try_variable;
- char ch, *rstinp, *atom_ext, *res_nam;
- Boolean mod = false;
- size_t len = LenChars( at->atom_ext, size_ID );
- memcpy( &try_variable, at, sizeof( Atom ) );
- atom_ext = try_variable.atom_ext;
- res_nam = try_variable.res_nam;
- rstinp = GetOptionCVal( opt_rstinp ); /* use value of restrinp option */
- /* OK if option menu has not been swapped */
- if ( EQn( rstinp, "xplor", 5 ) )
- {
- /* map XPLOR wildcard expression * to # */
- if ( atom_ext[len-1] == '*' )
- {
- atom_ext[len-1] = '#';
- mod = true;
- }
- /* map XPLOR wildcard expression % to # */
- else if ( atom_ext[len-1] == '%' )
- {
- atom_ext[len-1] = '#';
- mod = true;
- }
- /* map XPLOR wildcard expression + to # */
- else if ( atom_ext[len-1] == '+' )
- {
- atom_ext[len-1] = '#';
- mod = true;
- }
- /* map XPLOR wildcard expression # to ' ' and ## to # */
- else if ( atom_ext[len-1] == '#' )
- {
- atom_ext[len-1] = ' ';
- len--;
- mod = true;
- }
- }
- else if ( EQ( rstinp, "biosym" ) )
- {
- /* map BIOSYM wildcard expression * to X */
- if ( atom_ext[len-1] == '*' )
- {
- atom_ext[len-1] = 'X';
- mod = true;
- }
- /* map BIOSYM atom ...R to ...1 ... does not help if structure has 1HB etc. */
- /***disable
- else if ( atom_ext[len-1] == 'R' )
- {
- atom_ext[len-1] = '1';
- mod = true;
- }
- ***/
- /* map BIOSYM atom ...S to ...2 */
- /***disable
- else if ( atom_ext[len-1] == 'S' )
- {
- atom_ext[len-1] = '2';
- mod = true;
- }
- ***/
- /* corresponding maps for VAL HGR/S and LEU HDR/S */
- /***disable
- if ( ( memcmp( atom_ext, "HG", 2 ) == 0 &&
- memcmp( res_nam, "VAL", 3 ) == 0 ) ||
- ( memcmp( atom_ext, "HD", 2 ) == 0 &&
- memcmp( res_nam, "LEU", 3 ) == 0 ) )
- {
- if ( atom_ext[2] == 'R' )
- {
- atom_ext[2] = '1';
- mod = true;
- }
- else if ( atom_ext[2] == 'S' )
- {
- atom_ext[2] = '2';
- mod = true;
- }
- }
- ***/
- }
- /* map ILE HD to HD1 */
- if ( memcmp( atom_ext, "HD", 2 ) == 0 &&
- memcmp( res_nam, "ILE", 3 ) == 0 )
- {
- ch = atom_ext[2];
- if ( ch == ' ' ) /* because # may have been removed above */
- ch = '#';
- if ( ch != '1' )
- {
- memcpy( atom_ext, "HD1", 3 );
- atom_ext[3] = ch;
- mod = true;
- goto end;
- }
- }
- /* map ASN HD to HD2 */
- else if ( memcmp( atom_ext, "HD", 2 ) == 0 &&
- memcmp( res_nam, "ASN", 3 ) == 0 )
- {
- ch = atom_ext[2];
- if ( ch == ' ' ) /* because # may have been removed above */
- ch = '#';
- if ( ch != '2' )
- {
- memcpy( atom_ext, "HD2", 3 );
- atom_ext[3] = ch;
- mod = true;
- goto end;
- }
- }
- /* map GLN HE to HE2 */
- else if ( memcmp( atom_ext, "HE", 2 ) == 0 &&
- memcmp( res_nam, "GLN", 3 ) == 0 )
- {
- ch = atom_ext[2];
- if ( ch == ' ' ) /* because # may have been removed above */
- ch = '#';
- if ( ch != '2' )
- {
- memcpy( atom_ext, "HE2", 3 );
- atom_ext[3] = ch;
- mod = true;
- goto end;
- }
- }
- if ( mod ) goto end;
- /* map H to HN */
- if ( memcmp( atom_ext, "H ", 2 ) == 0 )
- {
- memcpy( atom_ext, "HN", 2 );
- mod = true;
- }
- /* map HN to H */
- else if ( memcmp( atom_ext, "HN", 2 ) == 0 )
- {
- memcpy( atom_ext, "H ", 2 );
- mod = true;
- }
- /* enable H..1/2/3 mapping */
- else if ( atom_ext[0] == 'H' )
- {
- ch = atom_ext[len-1];
- if ( ch == '1' || ch == '2' || ch == '3' )
- {
- atom_ext[0] = ch;
- memcpy( &(atom_ext[1]), at->atom_ext, len-1 );
- mod = true;
- }
- }
- else if ( atom_ext[1] == 'H' )
- {
- ch = atom_ext[0];
- if ( ch == '1' || ch == '2' || ch == '3' )
- {
- memcpy( atom_ext, &(at->atom_ext[1]), len-1 );
- atom_ext[len-1] = ch;
- mod = true;
- }
- }
- end:
- if ( ! mod )
- return( NULL );
- else
- return( &try_variable );
- }
- /**************** CmpResidOfAtomptr ****************/
- /*
- **
- */
- static int CmpResidOfAtomptr( const void *a, const void *b )
- {
- return( strcmp( (*(Atomptr *)a)->res_id , (*(Atomptr *)b)->res_id ) );
- }
- /**************** ExpandAtoms ****************/
- /*
- ** Expand the Atom list from 'wildcard' residue numbers to individual
- ** residues.
- **
- ** NOTE: res_num is treated as res_num_orig
- ** chain id not yet handled
- ** Handling is such that for ** wildcarded residues **
- ** occurrences in every chain are handled, and the residue info
- ** (including ID) is completed ** by copying from the residue list **.
- */
- static Mol *ExpandAtoms( Atom *atoms, const int at_count, Mol *molptr )
- {
- char *fct_name = "ExpandAtoms";
- int count, new_count = 0, res_count;
- Atom *list, *new_list, *new_atoms;
- /* Mol *new_molecule; */
- Residue *resptr;
- /* determine how many atoms are in the molecule for each 'wildcarded' residue number */
- /* use same loop mechanism as in actual setup of list (cf. below) */
- res_count = molptr->res_count;
- resptr = molptr->residues;
- while ( res_count-- )
- {
- list = atoms;
- count = at_count;
- while ( count-- )
- {
- if ( NOTEQID( resptr->res_nam, list->res_nam ) )
- {
- list++;
- continue;
- }
- if ( NOTEQID( resptr->chain, list->chain ) && list->chain[0] != '*' )
- {
- list++;
- continue;
- }
- if ( resptr->res_num_orig == list->res_num_orig ||
- list->res_num_orig == wild_RESNUM )
- {
- new_count++;
- }
- list++;
- }
- resptr++;
- }
- /* now setup the new list */
- /* Loop through the residues of the molptr molecule;
- ** for each residue loop through ALL atoms in the atoms list,
- ** and copy the atoms for which the residue name matches.
- ** The new list has the same residue order as the molptr,
- ** and within each residue the same atom order as the original list.
- */
- new_atoms = CREATE_NEW( Atom, new_count );
- if ( !new_atoms )
- NOMEM( "new_atoms" );
- new_list = new_atoms;
- res_count = molptr->res_count;
- resptr = molptr->residues;
- while ( res_count-- )
- {
- list = atoms;
- count = at_count;
- while ( count-- )
- {
- if ( NOTEQID( resptr->res_nam, list->res_nam ) )
- {
- list++;
- continue;
- }
- if ( NOTEQID( resptr->chain, list->chain ) && list->chain[0] != '*' )
- {
- list++;
- continue;
- }
- if ( resptr->res_num_orig == list->res_num_orig ||
- list->res_num_orig == wild_RESNUM )
- /* specific defs MUST come BEFORE wildcards in order to be found first by FindAtom */
- {
- /* copy the atom description */
- memcpy( new_list, list, sizeof( Atom ) );
- /* fill in the residue info */
- memcpy( new_list->chain, resptr->chain, size_ID );
- strcpy( new_list->res_id, resptr->res_id );
- new_list->res_insert = resptr->res_insert;
- new_list->res_num_orig = resptr->res_num_orig;
- new_list++;
- }
- list++;
- }
- resptr++;
- }
- /* WriteAtoms( stdout, new_atoms, new_count ); */
- /* return the new list */
- return( MakeDefsMol( new_atoms, new_count, true ) );
- }
- /**************** MakeDefsMol ****************/
- /*
- ** Create a dummy set and molecule, and make the residues for the new atom list
- ** (residues needed for use in FindAtom2)
- */
- static Mol *MakeDefsMol( Atom *new_atoms, const int new_count, const Boolean makeres )
- {
- char *fct_name = "MakeDefsMol";
- Strucset set; /* this is not a proper Strucset */
- Mol *new_molecule;
- new_molecule = CREATE_NEW( Mol, 1 );
- if ( !new_molecule )
- NOMEM( "new_molecule" );
- strcpy( set.name, "nameDefs" ); /* only these two needed... */
- set.molecule = new_molecule; /* ...in MakeResidues */
- new_molecule->atoms = new_atoms;
- new_molecule->atom_count = new_count;
-
- if ( makeres )
- MakeResidues( &set, false ); /* keep the original ordering! */
- /*
- ** Residue order already as in molptr, atom order same as in original list
- ** FindAtom2 executes "map on first match", therefore order is important
- ** Only if residues are scrambled FindAtom2 will fail
- */
- return( new_molecule );
- }