##################################
#                                #
# Last modified 01/24/2013       # 
#                                #
# Georgi Marinov                 #
#                                # 
##################################

import sys
import string
import scipy.stats
from sets import Set
import time

def run():

    if len(sys.argv) < 4:
        print 'usage: python %s list_of_files min_counts/copies p_value_cutoff outfile' % sys.argv[0]
        print '\tlist_of_files format:'
        print '\t#label\tfile_name\tgeneID_field\tp_value_fields\tcounts/copies fields\tfractionFieldID'
        print '\tyou can specify multiple p-value fields, comma-separated; all of them will have to pass the p-value threshold'
        print '\tthe counts/copies fields should include both maternal and paternal reads/copies'
        print '\tif you want to threshold on both read and copies, use the following separation: reads1,reads2|reads'
        sys.exit(1)

    list_of_files=sys.argv[1]
    minCopiesOrCounts = float(sys.argv[2])
    p_value_cutoff = float(sys.argv[3])
    outfilename = sys.argv[4]

    GeneDict={}

    Samples = []

    linelist=open(list_of_files)
    for line1 in linelist:
        fields1=line1.strip().split('\t')
        label = fields1[0]
        Samples.append(label)
        file = fields1[1]
        geneIDField = int(fields1[2])
        pValueIDFields = []
        fields = fields1[3].split(',')
        for ID in fields:
            pValueIDFields.append(int(ID))
        pValueIDFields.sort()
        CopiesOrCountsIDFields = []
        fields2 = fields1[4].split('|')
        for ID in fields2:
            countsfields = ID.split(',')
            fields_tuple = []
            for ID2 in countsfields:
                fields_tuple.append(int(ID2))
            CopiesOrCountsIDFields.append(tuple(fields_tuple))
        print label, CopiesOrCountsIDFields
        CopiesOrCountsIDFields .sort()
        fractionFieldID = int(fields1[5])
        lines = open(file)
        for line in lines:
            if line.startswith('#') or line.startswith('tracking_id'):
                continue
            fields=line.strip().split('\t')
            gene = fields[geneIDField]
            DoesNotPassCountsAndCopies = False
            for fields_tuple in CopiesOrCountsIDFields:
                CopiesOrCounts = 0
                for ID in list(fields_tuple):
                    CopiesOrCounts += float(fields[ID])
                if CopiesOrCounts < minCopiesOrCounts:
                    DoesNotPassCountsAndCopies = True
            if DoesNotPassCountsAndCopies:
                continue
            pvalues = []
            for ID in pValueIDFields:
                if fields[ID] == 'not_calculated':
                    pvalues.append(1)
                else:
                    pvalues.append(float(fields[ID]))
            if max(pvalues) > p_value_cutoff:
                continue
            if GeneDict.has_key(gene):
                pass
            else:
                GeneDict[gene]={}
            GeneDict[gene][label]=float(fields[fractionFieldID])

    outfile = open(outfilename, 'w')

    Samples.sort()

    outline = '#'
    for label in Samples:
        outline = outline + '\t' + label
    outfile.write(outline + '\n')

    genes = GeneDict.keys()
    genes.sort()
    for gene in genes:
        outline = gene
        for label in Samples:
            if GeneDict[gene].has_key(label):
                outline = outline + '\t' + str(GeneDict[gene][label])
            else:
                outline = outline + '\tnan'
        outfile.write(outline + '\n')

    outfile.close()

run()