+ moved penalty lookup and compute_doncc to Plif module
[qpalma.git] / qpalma / Plif.py
index a52594f..9fe0598 100644 (file)
@@ -1,8 +1,14 @@
 #!/usr/bin/env python
 # -*- coding: utf-8 -*-
 
+import math
+import pdb
+
+from numpy.matlib import zeros,isinf
+
 import QPalmaDP
 
+
 class Plf: 
    """
    means piecewise linear function
@@ -82,3 +88,66 @@ def logspace(a,b,n):
 
 def log10(x):
    return math.log(x)/math.log(10)
+
+
+def penalty_lookup_new(penalty_struct, value):
+
+   limits = penalty_struct.limits
+   penalties = penalty_struct.penalties
+
+   if penalty_struct.transform == 'log':
+      value = math.log(value)
+
+   elif penalty_struct.transform == 'log(+3)':
+      value = math.log(value+3)
+
+   elif penalty_struct.transform == 'log(+1)':
+      value = math.log(value+1)
+
+   elif penalty_struct.transform == '':
+      pass
+
+   elif penalty_struct.transform == '(+3)':
+      value += 3
+         
+   elif penalty_struct.transform == 'mod3':
+     value = value%3
+
+   else:
+      assert False,'unknown penalty transform'
+
+   count = len([elem for elem in limits if elem <= value])
+
+   if count == 0:
+      pen = penalties[0]
+   elif count == len(limits):
+      pen=penalties[-1]
+   else:
+      pen = (penalties[count]*(value-limits[count-1]) + penalties[count-1]\
+         *(limits[count]-value)) / (limits[count]-limits[count-1])
+
+   return pen
+
+
+def compute_donacc(donor_supp, acceptor_supp, d, a):
+
+   assert(len(donor_supp)==len(acceptor_supp))
+   size = len(donor_supp)
+  
+   donor    = [0.0]*size
+   acceptor = [0.0]*size
+
+   for idx in range(size):
+     if isinf(donor_supp[idx]):
+       donor[idx] = donor_supp[idx]
+     else:
+       donor[idx] = penalty_lookup_new(d, donor_supp[idx])
+
+     if isinf(acceptor_supp[idx]):
+       acceptor[idx] = acceptor_supp[idx]
+     else:
+       acceptor[idx] = penalty_lookup_new(a,acceptor_supp[idx])
+
+   return donor,acceptor
+
+