diff --git a/include/kitty/threshold_identification.hpp b/include/kitty/threshold_identification.hpp
index 02533d8e..d34953e2 100755
--- a/include/kitty/threshold_identification.hpp
+++ b/include/kitty/threshold_identification.hpp
@@ -30,47 +30,173 @@
   \author CS-472 2020 Fall students
 */
 
+
 #pragma once
 
+#include <stdio.h>
+#include <stdlib.h>
+#include <iostream>
+#include <cstdlib>
+#include <string>
+#include <fstream>
+#include <iostream>
+#include <cassert>
 #include <vector>
-// #include <lpsolve/lp_lib.h> /* uncomment this line to include lp_solve */
+#include <lpsolve/lp_lib.h> /* uncomment this line to include lp_solve */
+#include "isop.hpp"
 #include "traits.hpp"
+#include "constructors.hpp"
+#include "operations.hpp"
+#include "properties.hpp"
+#include "cube.hpp"
+#include "implicant.hpp"
+#include "bit_operations.hpp"
+#include "print.hpp"
 
-namespace kitty
-{
 
-/*! \brief Threshold logic function identification
 
-  Given a truth table, this function determines whether it is a threshold logic function (TF)
-  and finds a linear form if it is. A Boolean function is a TF if it can be expressed as
 
-  f(x_1, ..., x_n) = \sum_{i=1}^n w_i x_i >= T
 
-  where w_i are the weight values and T is the threshold value.
-  The linear form of a TF is the vector [w_1, ..., w_n; T].
+namespace kitty {
+
+    template<typename TT, typename = std::enable_if_t <kitty::is_complete_truth_table<TT>::value>>
+    bool is_threshold(const TT &tt, std::vector <int64_t> *plf = nullptr) {
+
+        std::vector <int64_t> linear_form;
+        std::vector <int8_t> negative_unate; 
+        TT tt0 = tt;
+        for (auto i = 0u; i < tt.num_vars(); i++) {
+
+            bool neg_unate = is_negative_unate_in_i(tt, i);
+            if (is_positive_unate_in_i(tt, i) + neg_unate == 0) {
+                return false;
+            }
+            if (neg_unate) {
+                flip_inplace(tt0, i);
+                negative_unate.push_back(i);
+            }
+
+        }
+        std::vector <kitty::cube> positiveConst = isop(tt0);
+        std::vector <kitty::cube> negativeConst = isop(unary_not(tt0));
+
+        lprec *lp;
+        int Numcol, *colnb = NULL, ret = 0;
+        REAL *row = NULL;
+        Numcol = tt.num_vars() + 1;
+
+        lp = make_lp(0, Numcol);
+        if (lp == NULL) {
+            return false; 
+        }
+
+        colnb = (int *) malloc(Numcol * sizeof(*colnb));
+        row = (REAL *) malloc(Numcol * sizeof(*row));
+
+        if ((colnb == NULL) || (row == NULL)) {
+            return false;
+        }
+        if (ret == 0) {
+            for (unsigned int j = 0; j <= tt.num_vars(); j++) {
+                colnb[0] = j + 1;
+                row[j] = 1;
+                add_constraintex(lp, 1, row, colnb, GE, 0);
+            }
+        }
+        if (ret == 0) {
+            for (auto k = 0; k < positiveConst.size(); k++) {
+                kitty::cube pos = positiveConst.at(k);
+                for (auto i = 0; i < tt.num_vars(); i++) {
+                    bool cond = pos.get_mask(i) && pos.get_bit(i); 
+                    colnb[i] = i + 1;
+                    if (cond == 1) {
+                        row[i] = 1;
+                    } else {
+                        row[i] = 0;
+                    }
+                }
+                colnb[tt.num_vars()] = tt.num_vars() + 1; 
+                row[tt.num_vars()] = -1.0;
+                add_constraintex(lp, Numcol, row, colnb, GE, 0);
+            }
+        }
+        if (ret == 0) {
+            for (auto k = 0; k < negativeConst.size(); k++) {
+                kitty::cube neg = negativeConst.at(k);
+                for (auto i = 0; i < tt.num_vars(); i++) {
+                    colnb[i] = i + 1;
+                    bool cond = !neg.get_mask(i) || (neg.get_mask(i) && neg.get_bit(i)); 
+                    if (cond) {
+                        row[i] = 1.0;
+                    } else {
+                        row[i] = 0.0;
+                    }
+                }
+                colnb[tt.num_vars()] = tt.num_vars() + 1;
+                row[tt.num_vars()] = -1.0;
+                add_constraintex(lp, Numcol, row, colnb, LE, -1);
+            }
+        }
+        if (ret == 0) {
+            set_add_rowmode(lp, FALSE); 
+            for (int j = 0; j < Numcol; j++) {
+                colnb[j] = j + 1;
+                row[j] = 1;
+            }
+            set_obj_fnex(lp, Numcol, row, colnb);
+        }
+        set_minim(lp);
+        set_verbose(lp, IMPORTANT);
+        ret = solve(lp);
+        if (ret != 0) {
+            return false;
+        } else {
+            get_variables(lp, row);
+            for (int i = 0; i < Numcol; i++) {
+                linear_form.push_back((int64_t) row[i]);
+            }
+            for (auto i : negative_unate) {
+                linear_form.at(i) = -linear_form.at(i);
+                linear_form[Numcol - 1] = linear_form[Numcol - 1] + linear_form[i];
+            }
+            *plf = linear_form;
+            if (colnb != NULL)
+                free(colnb);
+            if (row != NULL)
+                free(row);
+            if (lp != NULL)
+                delete_lp(lp);
+            return true;
+        }
+    }
+    
+    template<typename TT, typename = std::enable_if_t <kitty::is_complete_truth_table<TT>::value>>
+    bool is_negative_unate_in_i(const TT &tt, uint8_t i) {
+        auto const tt0 = cofactor0(tt, i);
+        auto const tt1 = cofactor1(tt, i);
+        for (auto bit = 0; bit < (2 << (tt.num_vars() - 1)); bit++) {
+            if (get_bit(tt0, bit) >= get_bit(tt1, bit)) {
+                continue;
+            } else {
+                return false;
+            }
+        }
+        return true;
+    }
+
+    template<typename TT, typename = std::enable_if_t <kitty::is_complete_truth_table<TT>::value>>
+    bool is_positive_unate_in_i(const TT &tt, uint8_t i) {
+        auto const tt0 = cofactor0(tt, i);
+        auto const tt1 = cofactor1(tt, i);
+        for (auto bit = 0; bit < (2 << (tt.num_vars() - 1)); bit++) {
+            if (get_bit(tt0, bit) <= get_bit(tt1, bit)) {
+                continue;
+            } else {
+                return false;
+            }
+        }
+        return true;
+    }
 
-  \param tt The truth table
-  \param plf Pointer to a vector that will hold a linear form of `tt` if it is a TF.
-             The linear form has `tt.num_vars()` weight values and the threshold value
-             in the end.
-  \return `true` if `tt` is a TF; `false` if `tt` is a non-TF.
-*/
-template<typename TT, typename = std::enable_if_t<is_complete_truth_table<TT>::value>>
-bool is_threshold( const TT& tt, std::vector<int64_t>* plf = nullptr )
-{
-  std::vector<int64_t> linear_form;
-
-  /* TODO */
-  /* if tt is non-TF: */
-  return false;
-
-  /* if tt is TF: */
-  /* push the weight and threshold values into `linear_form` */
-  if ( plf )
-  {
-    *plf = linear_form;
-  }
-  return true;
 }
 
-} /* namespace kitty */