Adding a contact file feature

src/INTBuilder.c

@@ -95,6 +95,22 @@ void getCandidatesForP1(struct pdb_values* pdb2, struct residue** t_candid1, str
     }
 }
 
+
+int check_ledger(int* ledger, int ledger_size, int entry) {
+    /* Checks the ledger array for the entry value.
+     * Returns 1 if present, else 0
+     */
+    int i = 0;
+
+    for (i = 0; i < ledger_size; i++) {
+        if(ledger[i] == entry){
+            return 1;
+        }
+    }
+
+    return 0;
+}
+
 void getInterface(struct pdb_values* pdbR, struct pdb_values* pdbL, int confID){
     /* This function will output the docking interface between pdbR and pdbL
      * To compute the docking interface, we considered that any residue of P1 at or closer than 10A
@@ -122,6 +138,28 @@ void getInterface(struct pdb_values* pdbR, struct pdb_values* pdbL, int confID){
     struct residue** t_candidateR = NULL;
     int nbCandR = pdbR->nbRes, nbCandL = pdbL->nbRes;
 
+    /* Contacts are stored in the following way for each array.
+     * First is the ID of the particle for the receptor protein, then 
+     * one or more ID for the particles in contact from the ligand protein.
+     * We go to the next ID of the receptor protein when
+     * we encounter the -1 value
+     *
+     * In the case of residue resolution, only residue ID is added. 
+     * In atom resolution, two values are added every time, first residue 
+     * then atom.
+     */
+    int* pdbR_contacts_id = NULL;
+    int allocate_contacts_cpt = pdbR->nbAtom, contact_iter = 0, new_resi_check = 1;
+
+    pdbR_contacts_id = calloc(allocate_contacts_cpt, sizeof(int));
+    if (pdbR_contacts_id == NULL) {
+        perror("calloc");
+        exit(EXIT_FAILURE);
+    }
+
+     /* Keeps in memory for the current resi the ones already in contact on the lig */
+    int t_res_contact_ledger[pdbL->nbRes], ledger_size = 0;
+
     /** Here we reset the number of clashes as this is a new interface computation **/
     clash = 0;
 
@@ -199,11 +237,28 @@ void getInterface(struct pdb_values* pdbR, struct pdb_values* pdbL, int confID){
     /* For each candidate residue of the receptor */
     for(i=0;i<nbCandR;i++){
         isClashed = 0;
+        if (!atom_res && !new_resi_check){
+            new_resi_check = 1;
+            pdbR_contacts_id[contact_iter] = -1;
+            contact_iter++;
+            ledger_size = 0;
+
+            check_size_contacts(&pdbR_contacts_id, &allocate_contacts_cpt, contact_iter);
+        }
+
         for(j=0;j<t_candidateR[i]->nbAtom;j++){
             x1 = t_candidateR[i]->x[j];
             y1 = t_candidateR[i]->y[j];
             z1 = t_candidateR[i]->z[j];
 
+            if (atom_res && !new_resi_check){
+                new_resi_check = 1;
+                pdbR_contacts_id[contact_iter] = -1;
+                contact_iter++;
+
+                check_size_contacts(&pdbR_contacts_id, &allocate_contacts_cpt, contact_iter);
+            }
+
             /* For each candidate residue of the ligand */
             for(k=0;k<nbCandL;k++){
                 for(l=0;l<t_candidateL[k]->nbAtom;l++){
@@ -228,8 +283,38 @@ void getInterface(struct pdb_values* pdbR, struct pdb_values* pdbL, int confID){
                         t_candidateR[i]->isAtomCandidate[j] = 1;
                         t_candidateL[k]->isAtomCandidate[l] = 1;
 
-                        if(!atom_res)
+                        /* Check for contacts */
+                        if(atom_res){
+                            if(new_resi_check){
+                                pdbR_contacts_id[contact_iter] = i;
+                                pdbR_contacts_id[contact_iter+1] = j;
+                                contact_iter += 2;
+                                new_resi_check = 0;
+                                check_size_contacts(&pdbR_contacts_id, &allocate_contacts_cpt, contact_iter);
+                            }
+                            pdbR_contacts_id[contact_iter] = k;
+                            pdbR_contacts_id[contact_iter+1] = l;
+                            contact_iter += 2;
+                            check_size_contacts(&pdbR_contacts_id, &allocate_contacts_cpt, contact_iter);
+                        }else{
+                            if(new_resi_check){
+                                pdbR_contacts_id[contact_iter] = i;
+                                contact_iter++;
+                                new_resi_check = 0;
+
+                                check_size_contacts(&pdbR_contacts_id, &allocate_contacts_cpt, contact_iter);
+                            }
+
+                            /* Check if this contact has already been tagged for this residue */
+                            if(!check_ledger(t_res_contact_ledger, ledger_size, k)){
+                                pdbR_contacts_id[contact_iter] = k;
+                                t_res_contact_ledger[ledger_size] = k;
+                                ledger_size++;
+                                contact_iter++;
+                                check_size_contacts(&pdbR_contacts_id, &allocate_contacts_cpt, contact_iter);
+                            }
                             break;
+                        }
                     }
                 }
                 if(clash >= TOO_MANY_CLASHES){
@@ -246,6 +331,13 @@ void getInterface(struct pdb_values* pdbR, struct pdb_values* pdbL, int confID){
     }
 
 
+    if(pdbR_contacts_id[contact_iter-1] == -1){
+        pdbR_contacts_id[contact_iter-1] = -2;
+    }else{
+        pdbR_contacts_id[contact_iter] = -2;
+    }
+    contact_iter++;
+
     /* If the amount of clashes is still ok */
     if(clash < TOO_MANY_CLASHES){
         /* Output the ligand's interface */
@@ -257,10 +349,16 @@ void getInterface(struct pdb_values* pdbR, struct pdb_values* pdbL, int confID){
         fprintf(outputFile_rec,"%d ",confID);
         write_candidate(t_candidateR,nbCandR,outputFile_rec);
         fprintf(outputFile_rec,"\n");
+
+        /* Output the receptor's interface */
+        fprintf(outputFile_contacts,"%d ",confID);
+        write_contacts(t_candidateR,t_candidateL,pdbR_contacts_id,outputFile_contacts);
+        fprintf(outputFile_contacts,"\n");
     }
 
     free(t_candidateL);
     free(t_candidateR);
+    free(pdbR_contacts_id);
 }
 
 int main(int argc, char** argv){

src/INTBuilder.h

@@ -3,6 +3,7 @@
 #ifndef INTB_HEADER
 #define INTB_HEADER
 
+int check_ledger(int* ledger, int ledger_size, int entry);
 void getCandidatesForP1(struct pdb_values* pdb2, struct residue** t_candid1, struct residue** t_candid2, int* nbCand1, int* nbCand2);
 void getInterface(struct pdb_values* pdbR, struct pdb_values* pdbL, int confID);
 void removeSpace(char* string);

src/allocate.c

@@ -304,6 +304,19 @@ void free_argLine(){
     }
 }
 
+
+void check_size_contacts(int** t_contacts, int* allocated_size, int nb_contacts){
+    /* Checks if the allocated_size for t_contacts is enough for nb_contacts
+     * If not, doubles the allocated size
+     */
+
+    if (nb_contacts >= *allocated_size - 10) {
+        *allocated_size *= 2;
+        *t_contacts = realloc(*t_contacts, (*allocated_size)*sizeof(int));
+    }
+
+}
+
 void freeAll(struct pdb_values* pdbR, struct pdb_values* pdbL, struct docking_results* dock_res){
     /* Free all the arguments passed to this function */
     freePDB(pdbL);

src/allocate.h

@@ -10,5 +10,6 @@ void reset_dockingResults(struct docking_results* dock_res);
 struct pdb_values* allocate_pdb(int nbRes);
 void free_argLine();
 void freeAll(struct pdb_values* pdbR, struct pdb_values* pdbL, struct docking_results* dock_res);
+void check_size_contacts(int** t_contacts, int* allocated_size, int nb_contacts);
 
 #endif

src/fileIO.c

@@ -857,25 +857,10 @@ void write_candidate(struct residue** t_candid, int sizeCand, FILE* output_strea
                 if(t_candid[i]->isAtomCandidate[j]){
 
                     /* Getting the residue name */
-                    len = strlen(t_candid[i]->idRes);
-                    k = 0;
-                    for(l=0;l<len;l++){
-                        if(t_candid[i]->idRes[l] != ' '){
-                            buf[k] = t_candid[i]->idRes[l];
-                            k++;
-                        }
-                    }
-                    buf[k] = '\0';
+                    get_resi_name(t_candid, i, buf);
 
                     /* Now, getting the atom name */
-                    len = strlen(t_candid[i]->atomType[j]);
-                    k = 0;
-                    for(l=0;l<len;l++){
-                        if(t_candid[i]->atomType[j][l] != ' '){
-                            buf2[k] = t_candid[i]->atomType[j][l];
-                            k++;
-                    }            }
-                    buf2[k] = '\0';
+                    get_atom_name(t_candid, i, j, buf2);
                     fprintf(output_stream,"'%s' '%c' '%s' ",buf,t_candid[i]->chain, buf2);
                     fflush(output_stream);
                 }
@@ -899,3 +884,108 @@ void write_candidate(struct residue** t_candid, int sizeCand, FILE* output_strea
     }
 }
 
+
+void get_resi_name(struct residue** t_residues, int resi_id, char* str_resi){
+    int len = 0, k = 0, l = 0;
+
+    /* Getting the residue name */
+    len = strlen(t_residues[resi_id]->idRes);
+    k = 0;
+    for(l=0;l<len;l++){
+        if(t_residues[resi_id]->idRes[l] != ' '){
+            str_resi[k] = t_residues[resi_id]->idRes[l];
+            k++;
+        }
+    }
+    str_resi[k] = '\0';
+}
+
+
+void get_atom_name(struct residue** t_residues, int resi_id, int atom_id, char* str_atom){
+    int len = 0, k = 0, l = 0;
+
+    /* Now, getting the atom name */
+    len = strlen(t_residues[resi_id]->atomType[atom_id]);
+    k = 0;
+    for(l=0;l<len;l++){
+        if(t_residues[resi_id]->atomType[atom_id][l] != ' '){
+            str_atom[k] = t_residues[resi_id]->atomType[atom_id][l];
+            k++;
+    }            }
+    str_atom[k] = '\0';
+}
+
+
+void write_contacts(struct residue** t_candidR, struct residue** t_candidL, int* t_contacts_id, FILE* output_stream){
+    /* Output the residues selected as being part of the interface in the array
+     * t_candid
+     */
+    int i = 0, j = 0, k= 0, l = 0, len = 0;
+    int resi_id = 0, atom_id = 0, new_contact = 1;
+    char buf[15] = "", buf2[15] = "";
+
+    if(atom_res){
+        while(t_contacts_id[i] != -2) {
+            resi_id = t_contacts_id[i];
+            atom_id = t_contacts_id[i+1];
+
+            /* Value is -1 if we look at a new contact group */
+            if(resi_id == -1){
+                new_contact = 1;
+                i++;
+
+                resi_id = t_contacts_id[i];
+                atom_id = t_contacts_id[i+1];
+
+                fprintf(output_stream,"-- ");
+                fflush(output_stream);
+            }
+
+
+            /* Getting the residue and atom name */
+            if(new_contact){
+                new_contact = 0;
+                get_resi_name(t_candidR, resi_id, buf);
+                get_atom_name(t_candidR, resi_id, atom_id, buf2);
+                fprintf(output_stream,"'%s' '%c' '%s' ",buf,t_candidR[resi_id]->chain, buf2);
+            }else{
+                get_resi_name(t_candidL, resi_id, buf);
+                get_atom_name(t_candidL, resi_id, atom_id, buf2);
+                fprintf(output_stream,"'%s' '%c' '%s' ",buf,t_candidL[resi_id]->chain, buf2);
+            }
+            fflush(output_stream);
+
+            /* Check next contact */
+            i += 2;
+        }
+    }else{
+        while(t_contacts_id[i] != -2) {
+            resi_id = t_contacts_id[i];
+
+            /* Value is -1 if we look at a new contact group */
+            if(resi_id == -1){
+                new_contact = 1;
+                i++;
+
+                resi_id = t_contacts_id[i];
+
+                fprintf(output_stream,"-- ");
+                fflush(output_stream);
+            }
+
+            /* Getting the residue and atom name */
+            if(new_contact){
+                new_contact = 0;
+                get_resi_name(t_candidR, resi_id, buf);
+                fprintf(output_stream,"'%s' '%c' ",buf,t_candidR[resi_id]->chain);
+            }else{
+                get_resi_name(t_candidL, resi_id, buf);
+                fprintf(output_stream,"'%s' '%c' ",buf,t_candidL[resi_id]->chain);
+            }
+            fflush(output_stream);
+
+            /* Check next contact */
+            i++;
+        }
+    }
+}

src/fileIO.h

@@ -14,5 +14,8 @@ struct docking_results* getDataForComplex_HEX();
 struct docking_results* getDataForComplex_MAXDO();
 struct docking_results* getDataForComplex_HCMD2();
 struct docking_results* getDataForComplex_DEFAULT();
+void get_resi_name(struct residue** t_residues, int resi_id, char* str_resi);
+void get_atom_name(struct residue** t_residues, int resi_id, int atom_id, char* str_atom);
+void write_contacts(struct residue** t_candidR, struct residue** t_candidL, int* t_contacts_id, FILE* output_stream);
 
 #endif

src/param.c

@@ -277,6 +277,13 @@ struct argLine* parseLineOfArgument(int argc, char** argv){
                 perror("fopen");
                 exit(EXIT_FAILURE);
             }
+
+            sprintf(buf,"%s/%s-%s_rec_%d_contacts.txt",outputDir,receptor,ligand,target_conf);
+            outputFile_contacts = fopen(buf,"w");
+            if(outputFile_contacts == NULL){
+                perror("fopen");
+                exit(EXIT_FAILURE);
+            }
         }else if(complexPDB){
             sprintf(buf,"%s/%s-%s_rec_0_dockinter.txt",outputDir,receptor,ligand);
             outputFile_rec = fopen(buf,"w");
@@ -291,6 +298,13 @@ struct argLine* parseLineOfArgument(int argc, char** argv){
                 perror("fopen");
                 exit(EXIT_FAILURE);
             }
+
+            sprintf(buf,"%s/%s-%s_rec_0_contacts.txt",outputDir,receptor,ligand);
+            outputFile_contacts = fopen(buf,"w");
+            if(outputFile_contacts == NULL){
+                perror("fopen");
+                exit(EXIT_FAILURE);
+            }
         }else{
             sprintf(buf,"%s/%s-%s_rec_dockinter.txt",outputDir,receptor,ligand);
             outputFile_rec = fopen(buf,"w");
@@ -305,6 +319,13 @@ struct argLine* parseLineOfArgument(int argc, char** argv){
                 perror("fopen");
                 exit(EXIT_FAILURE);
             }
+
+            sprintf(buf,"%s/%s-%s_rec_contacts.txt",outputDir,receptor,ligand);
+            outputFile_contacts = fopen(buf,"w");
+            if(outputFile_contacts == NULL){
+                perror("fopen");
+                exit(EXIT_FAILURE);
+            }
         }
     }
 

src/struct.h

@@ -68,6 +68,7 @@ char*   outputDir;          /* Output directory for the interface files */
 char*   outputPDB;          /* Output directory if we construct the PDB(s) */
 FILE*   outputFile_rec;     /* Name of the interface output file for the receptor */
 FILE*   outputFile_lig;     /* Name of the interface output file for the ligand */
+FILE*   outputFile_contacts;/* Name of the output file for the contacts */
 
 /*************************************************/