From 91915e5fa8c897fffdccc7b74fed57d69e0d8489 Mon Sep 17 00:00:00 2001
From: tforest <thomas.forest@college-de-france.fr>
Date: Thu, 18 Nov 2021 11:14:50 +0100
Subject: [PATCH] Modifs script C

---
 vcf_to_sfs.c | 65 ++++++++++++++++++++++++++++++++++++++++++++++++++--
 1 file changed, 63 insertions(+), 2 deletions(-)

diff --git a/vcf_to_sfs.c b/vcf_to_sfs.c
index 8ea4c8d..3581e2e 100644
--- a/vcf_to_sfs.c
+++ b/vcf_to_sfs.c
@@ -19,8 +19,42 @@ void slice_str(const char * str, char * buffer, size_t start, size_t end)
     buffer[j] = 0;
 }
 
+int min(int * array, int size){
+    //Consider first element as smallest
+   int smallest = array[0];
+   int i;
+   for (i = 0; i < num; i++) {
+      if (a[i] < smallest) {
+         smallest = a[i];
+      }
+   }
+}
+
+int countDistinct(int a[], int n)      //Function Definition
+{
+   int i, j, count = 0;
+   //Traverse the array
+   for (i = 1; i < n; i++)      //hold an array element
+   {
+      for (j = 0; j < i; j++)   
+      {
+         if (a[i] == a[j])    //Check for duplicate elements
+         {
+            break;             //If duplicate elements found then break
+         }
+      }
+      if (i == j)
+      {
+         count++;     //increment the number of distinct elements
+      }
+   }
+   return count;      //Return the number of distinct elements
+}
 
 # define LL 8192   /* line length maximum */
+# define DIPLOID true
+# define FOLDED true
+# define IGNORED_FIELDS 9
 
 int main ( int argc, char *argv[] ){
     if ( argc < 3) {
@@ -29,9 +63,20 @@ int main ( int argc, char *argv[] ){
     }
     gzFile fp;
     char line[LL];
+    int N;
     char delim[] = "\t";
     fp = gzopen( argv[1], "r" );
 
+    // pop of size 2N when diploid
+    if (DIPLOID == true && FOLDED == false) {
+	N = 2 * atoi(argv[2]);
+	    } else {
+	N = atoi(argv[2]);
+    }
+
+    int snp_genotypes[N];
+    int SFS_values[N];
+
     gzgets( fp, line, LL );
     while ( ! gzeof( fp ) ){
 	int k = 0;
@@ -43,14 +88,30 @@ int main ( int argc, char *argv[] ){
 	char *vcf_field = strtok(line, delim);
 	while(vcf_field != NULL){
 	    k++;
-	    if (k > 9) {
+	    if (k > IGNORED_FIELDS) {
 		const size_t len = strlen(vcf_field);
 		char buffer[len + 1];
 		//printf("'%s'\n", ptr);
 		slice_str(vcf_field, buffer, 0, 0);
-		printf("%s      ", buffer);
+		//printf("%d %s      ", N, buffer);
+		snp_genotypes[k-IGNORED_FIELDS] = atoi(buffer);
+		//printf("%d ", smpl_genotype[k-9]);
 	    }
 	    vcf_field = strtok(NULL, delim);
+	    int c= countDistinct(snp_genotypes, N);
+            // skip if all individuals have the same genotype
+	    if (c == 1) {
+		continue;
+		gzgets( fp, line, LL );
+	    }
+	    /* int i; */
+	    /* for (i = 1; i < N; ++i) */
+	    /* 	{ */
+	    /* 	    printf("%d ", snp_genotypes[i]); */
+	    /* 	} */
+	    int allele_counts[c];
+	    
+	    min(allele_counts, N);
 	}
 	// printf("%s", line );
 	// loads the next line