# Path: nanopore/analyses/utils.py
# class UniqueList(collections.MutableSet):
# class AlignedPair:
#     def __init__(self, iterable=None):
#     def __len__(self):
#     def __contains__(self, key):
#     def add(self, key):
#     def discard(self, key):
#     def __iter__(self):
#     def __reversed__(self):
#     def pop(self, last=True):
#     def remove(self, last=True):
#     def __getitem__(self, index):
#     def __repr__(self):
#     def __eq__(self, other):
# def pathToBaseNanoporeDir():
#     def __init__(self, refPos, refSeq, readPos, isReversed, readSeq, pPair):
#     def isMatch(self):
#     def isMismatch(self):
#     def getRefBase(self):
#     def getReadBase(self):
#     def getSignedReadPos(self):
#     def getPrecedingReadInsertionLength(self, globalAlignment=False):
#     def getPrecedingReadDeletionLength(self, globalAlignment=False):
#     def _indelLength(pos, pPos):
#     def iterator(alignedRead, refSeq, readSeq):
# def getAbsoluteReadOffset(alignedRead, refSeq, readSeq):
# def getExonerateCigarFormatString(alignedRead, sam):
# def samToBamFile(samInputFile, bamOutputFile):
# def getFastaDictionary(fastaFile):
# def getFastqDictionary(fastqFile):
# def makeFastaSequenceNamesUnique(inputFastaFile, outputFastaFile):
# def makeFastqSequenceNamesUnique(inputFastqFile, outputFastqFile):
# def normaliseQualValues(inputFastqFile, outputFastqFile):
# def samIterator(sam):
# def mergeChainedAlignedReads(chainedAlignedReads, refSequence, readSequence):
# def chainFn(alignedReads, refSeq, readSeq, scoreFn=lambda alignedRead, refSeq, readSeq : len(list(AlignedPair.iterator(alignedRead, refSeq, readSeq))), maxGap=200):
#     def getStartAndEndCoordinates(alignedRead):
# def combineSamFiles(baseSamFile, extraSamFiles, outputSamFile):
# def chainSamFile(samFile, outputSamFile, readFastqFile, referenceFastaFile, chainFn=chainFn):
# def learnModelFromSamFileTargetFn(target, samFile, readFastqFile, referenceFastaFile, outputModel):
# def learnModelFromSamFileTargetFn2(target, unnormalisedOutputModel, outputModel):
# def realignSamFileTargetFn(target, samFile, outputSamFile, readFastqFile, 
#                            referenceFastaFile, gapGamma, matchGamma, hmmFile=None, trainHmmFile=False, chainFn=chainFn):
# def realignSamFile2TargetFn(target, samFile, outputSamFile, readFastqFile, referenceFastaFile, hmmFile, gapGamma, matchGamma):
# def realignCigarTargetFn(target, exonerateCigarString, referenceSequenceName, referenceSequence, querySequenceName, querySequence, outputCigarFile, hmmFile, gapGamma, matchGamma):
# def realignSamFile3TargetFn(target, samFile, outputSamFile, tempCigarFiles):
# def normaliseHmmByReferenceGCContent(hmm, gcContent):
# def modifyHmmEmissionsByExpectedVariationRate(hmm, substitutionRate):
# def setHmmIndelEmissionsToBeFlat(hmm):


from cactus.bar.cactus_expectationMaximisation import Hmm, SYMBOL_NUMBER
from nanopore.analyses.utils import setHmmIndelEmissionsToBeFlat, normaliseHmmByReferenceGCContent, modifyHmmEmissionsByExpectedVariationRate, toMatrix
import sys
import numpy as np
###Modify HMMs output from EM training to normalise for background nucleotide frequencies.

def main():
    print "ARGS", sys.argv
    #Load HMM
    hmm = Hmm.loadHmm(sys.argv[1])

    #setHmmIndelEmissionsToBeFlat(hmm)

    #Normalise background emission frequencies, if requested to GC% given
    gcContent = float(sys.argv[2])
    print "Got GC content", gcContent
    normaliseHmmByReferenceGCContent(hmm, gcContent)
    
    #Modify match emissions by proposed substitution rate
    substitutionRate = float(sys.argv[3])
    print "Got substitution rate", substitutionRate
    modifyHmmEmissionsByExpectedVariationRate(hmm, substitutionRate)
    
    for state in range(0, hmm.stateNumber):
        n = toMatrix(hmm.emissions[(SYMBOL_NUMBER**2) * state:(SYMBOL_NUMBER**2) * (state+1)])
        print "For state, ref frequencies", map(sum, n)
        print "For state, read frequencies", map(sum, np.transpose(n))
    
    #Write out HMM
    hmm.write(sys.argv[4])

if __name__ == '__main__':

# Path: nanopore/analyses/utils.py
# class UniqueList(collections.MutableSet):
# class AlignedPair:
#     def __init__(self, iterable=None):
#     def __len__(self):
#     def __contains__(self, key):
#     def add(self, key):
#     def discard(self, key):
#     def __iter__(self):
#     def __reversed__(self):
#     def pop(self, last=True):
#     def remove(self, last=True):
#     def __getitem__(self, index):
#     def __repr__(self):
#     def __eq__(self, other):
# def pathToBaseNanoporeDir():
#     def __init__(self, refPos, refSeq, readPos, isReversed, readSeq, pPair):
#     def isMatch(self):
#     def isMismatch(self):
#     def getRefBase(self):
#     def getReadBase(self):
#     def getSignedReadPos(self):
#     def getPrecedingReadInsertionLength(self, globalAlignment=False):
#     def getPrecedingReadDeletionLength(self, globalAlignment=False):
#     def _indelLength(pos, pPos):
#     def iterator(alignedRead, refSeq, readSeq):
# def getAbsoluteReadOffset(alignedRead, refSeq, readSeq):
# def getExonerateCigarFormatString(alignedRead, sam):
# def samToBamFile(samInputFile, bamOutputFile):
# def getFastaDictionary(fastaFile):
# def getFastqDictionary(fastqFile):
# def makeFastaSequenceNamesUnique(inputFastaFile, outputFastaFile):
# def makeFastqSequenceNamesUnique(inputFastqFile, outputFastqFile):
# def normaliseQualValues(inputFastqFile, outputFastqFile):
# def samIterator(sam):
# def mergeChainedAlignedReads(chainedAlignedReads, refSequence, readSequence):
# def chainFn(alignedReads, refSeq, readSeq, scoreFn=lambda alignedRead, refSeq, readSeq : len(list(AlignedPair.iterator(alignedRead, refSeq, readSeq))), maxGap=200):
#     def getStartAndEndCoordinates(alignedRead):
# def combineSamFiles(baseSamFile, extraSamFiles, outputSamFile):
# def chainSamFile(samFile, outputSamFile, readFastqFile, referenceFastaFile, chainFn=chainFn):
# def learnModelFromSamFileTargetFn(target, samFile, readFastqFile, referenceFastaFile, outputModel):
# def learnModelFromSamFileTargetFn2(target, unnormalisedOutputModel, outputModel):
# def realignSamFileTargetFn(target, samFile, outputSamFile, readFastqFile, 
#                            referenceFastaFile, gapGamma, matchGamma, hmmFile=None, trainHmmFile=False, chainFn=chainFn):
# def realignSamFile2TargetFn(target, samFile, outputSamFile, readFastqFile, referenceFastaFile, hmmFile, gapGamma, matchGamma):
# def realignCigarTargetFn(target, exonerateCigarString, referenceSequenceName, referenceSequence, querySequenceName, querySequence, outputCigarFile, hmmFile, gapGamma, matchGamma):
# def realignSamFile3TargetFn(target, samFile, outputSamFile, tempCigarFiles):
# def normaliseHmmByReferenceGCContent(hmm, gcContent):
# def modifyHmmEmissionsByExpectedVariationRate(hmm, substitutionRate):
# def setHmmIndelEmissionsToBeFlat(hmm):


# Path: scripts/modifyHmm.py
from cactus.bar.cactus_expectationMaximisation import Hmm, SYMBOL_NUMBER
from nanopore.analyses.utils import setHmmIndelEmissionsToBeFlat, normaliseHmmByReferenceGCContent, modifyHmmEmissionsByExpectedVariationRate, toMatrix
import sys
import numpy as np

###Modify HMMs output from EM training to normalise for background nucleotide frequencies.

def main():
    print "ARGS", sys.argv
    #Load HMM
    hmm = Hmm.loadHmm(sys.argv[1])

    #setHmmIndelEmissionsToBeFlat(hmm)

    #Normalise background emission frequencies, if requested to GC% given
    gcContent = float(sys.argv[2])
    print "Got GC content", gcContent
    normaliseHmmByReferenceGCContent(hmm, gcContent)
    
    #Modify match emissions by proposed substitution rate
    substitutionRate = float(sys.argv[3])
    print "Got substitution rate", substitutionRate
    modifyHmmEmissionsByExpectedVariationRate(hmm, substitutionRate)
    
    for state in range(0, hmm.stateNumber):
        n = toMatrix(hmm.emissions[(SYMBOL_NUMBER**2) * state:(SYMBOL_NUMBER**2) * (state+1)])
        print "For state, ref frequencies", map(sum, n)
        print "For state, read frequencies", map(sum, np.transpose(n))
    
    #Write out HMM
    hmm.write(sys.argv[4])

if __name__ == '__main__':

    main()