All SNPs computed

1 files changed, 33 insertions, 15 deletions

diff --git a/src/lmm.cpp b/src/lmm.cpp
index d015988..1319d7c 100644
--- a/src/lmm.cpp
+++ b/src/lmm.cpp
@@ -1868,7 +1868,7 @@ void LMM::Analyze(std::function< SnpNameValues(size_t) >& fetch_snp,
  */
 
 
-void LMM::mdb_analyze(std::function< SnpNameValues(size_t) >& fetch_snp,
+void LMM::mdb_analyze(std::function< SnpNameValues2(size_t) >& fetch_snp,
                       const gsl_matrix *U, const gsl_vector *eval,
                       const gsl_matrix *UtW, const gsl_vector *Uty,
                       const gsl_matrix *W, const gsl_vector *y,
@@ -1995,16 +1995,20 @@ void LMM::mdb_analyze(std::function< SnpNameValues(size_t) >& fetch_snp,
 
   assert(num_markers > 0);
   for (size_t t = 0; t < num_markers; ++t) {
+  // for (size_t t = 0; t < 2; ++t) {
     if (t % progress_step == 0 || t == (num_markers - 1)) {
       ProgressBar("Reading SNPs", t, num_markers - 1);
     }
     // if (indicator_snp[t] == 0)
     // continue;
 
-    auto tup = fetch_snp(t); // use the call back
-    auto snp = get<0>(tup);
-    auto gs = get<1>(tup);
-    cout << "SNP: " << snp << endl;
+    auto tup = fetch_snp(t); // use the callback
+    auto success = get<0>(tup);
+    if (!success)
+      continue;
+    auto snp = get<1>(tup);
+    auto gs = get<2>(tup);
+    cout << t << " SNP: " << snp << endl;
 
     // check whether SNP is included in gwasnps (used by LOCO)
     /*
@@ -2096,31 +2100,45 @@ void LMM::mdb_calc_gwa(const gsl_matrix *U, const gsl_vector *eval,
   auto num_markers = stat.ms_entries;
 
   // fetch_snp is a callback function for every SNP row
-  // returns typedef std::tuple<string,std::vector<double> > SnpNameValues;
+  // returns typedef std::tuple<bool,string,std::vector<double> > SnpNameValues;
 
   size_t prev_line = 0;
   auto mdb_fetch = MDB_FIRST;
 
   auto cursor = lmdb::cursor::open(rtxn, geno_mdb);
 
-  std::function<SnpNameValues(size_t)>  fetch_snp = [&](size_t num) {
+  std::function<SnpNameValues2(size_t)>  fetch_snp = [&](size_t num) {
 
     string_view key,value;
 
+    /*
     while (prev_line <= num) {
       cursor.get(key, value, MDB_NEXT);
       prev_line++;
     }
-    cursor.get(key, value, mdb_fetch);
+    */
+    auto result = cursor.get(key, value, mdb_fetch);
     mdb_fetch = MDB_NEXT;
-    size_t num_floats = value.size() / sizeof(float);
-    assert(num_floats == ni_total);
-    const float* gsbuf = reinterpret_cast<const float*>(value.data());
-    auto snp = string(key);
-    std::vector<double> gs(gsbuf,gsbuf+sizeof(float)*value.size());
-    cout << "!!!!" << snp << endl << endl;
 
-    return std::make_tuple(snp,gs);
+    string snp;
+    vector<double> gs;
+
+    if (result) {
+      size_t num_floats = value.size() / sizeof(float);
+      cout << string(key) << ":" << num_floats << "," << ni_total << " " << endl;
+      assert(num_floats == ni_total);
+      const float* gsbuf = reinterpret_cast<const float*>(value.data());
+      auto snp = string(key);
+      // std::vector<double> gs(gsbuf,gsbuf+sizeof(float)*value.size());
+      auto size = num_floats;
+      gs.reserve(size);
+      for (size_t i = 0; i < size; ++i) {
+        gs.push_back(static_cast<double>(gsbuf[i]));
+      }
+      cout << "!!!!" << size << snp << ": " << gs[0] << "," << gs[1] << "," << gs[2] << "," << gs[3] << endl;
+    }
+    // return { result,snp,gs };
+    return make_tuple(result, snp, gs);
   };
 
   LMM::mdb_analyze(fetch_snp,U,eval,UtW,Uty,W,y,gwasnps,num_markers);