merge develop branch

.github/workflows/make.yml

@@ -25,4 +25,7 @@ jobs:
       uses: actions/upload-artifact@v4
       with:
         name: test-output
-        path: test/output/
+        path: |
+          test/output
+          test/testAll.sh
+          test/testAll.sh.log

Makefile

@@ -1,144 +1,138 @@
-# make FLAGS=-D__REGRESSION__
-#c++ -I /opt/homebrew/Cellar/eigen/3.4.0_1/include/ -c -I/Users/fvr124/metaDMG-cpp/htslib   regression.cpp -std=c++14  -D__REGRESSION__
+# --- Flags og kompileringsvalg ---
+FLAGS     := -O2
+CFLAGS    := $(FLAGS)
+CXXFLAGS  := $(FLAGS)
+CPPFLAGS  := $(CPPFLAGS) -Wall -Wextra
+LDFLAGS   := -lgsl
 
 CC  ?= gcc
 CXX ?= g++
 
-LIBS = -lz -lm -lbz2 -llzma -lpthread -lcurl -lgsl -lgslcblas
-
-CRYPTO_TRY=$(shell echo 'int main(){}'|$(CXX) -x c++ - -lcrypto 2>/dev/null -o /dev/null; echo $$?)
-ifeq "$(CRYPTO_TRY)" "0"
-$(info Crypto library is available to link; adding -lcrypto to LIBS)
-LIBS += -lcrypto
-else
-$(info Crypto library is not available to link; will not use -lcrypto)
-endif
+CXXSRC := $(wildcard *.cpp)
+CSRC   := $(wildcard *.c)
+OBJ    := $(CSRC:.c=.o) $(CXXSRC:.cpp=.o)
 
 
-#if htslib source is defined
-ifdef HTSSRC
-
-#if hts source is set to systemwide
-ifeq ($(HTSSRC),systemwide)
-$(info HTSSRC set to systemwide; assuming systemwide installation)
-LIBS += -lhts
-
-else
-
-#if hts source path is given
-# Adjust $(HTSSRC) to point to your top-level htslib directory
-$(info HTSSRC defined: $(HTSSRC))
-CPPFLAGS += -I"$(realpath $(HTSSRC))"
-LIBHTS := $(HTSSRC)/libhts.a
-LIBS := $(LIBHTS) $(LIBS)
 
+# --- Brugervalgte paths og biblioteker ---
+ifdef PREFIX
+CPPFLAGS += -I$(PREFIX)/include
+LDFLAGS  := -L$(PREFIX)/lib $(LDFLAGS)
 endif
 
-#if htssrc not defined
-else
-
-$(info HTSSRC not defined; using htslib submodule)
-$(info Use `make HTSSRC=/path/to/htslib` to build metadamage using a local htslib installation)
-$(info Use `make HTSSRC=systemwide` to build metadamage using the systemwide htslib installation)
-
-
-HTSSRC := $(CURDIR)/htslib
-CPPFLAGS += -I$(HTSSRC)
-LIBHTS := $(HTSSRC)/libhts.a
-LIBS := $(LIBHTS) $(LIBS)
-
-all: .activate_module
-
+ifdef EXTRA_INC
+CPPFLAGS += $(addprefix -I,$(EXTRA_INC))
 endif
 
-.PHONY: .activate_module
-
-.activate_module:
-	git submodule update --init --recursive
-	$(MAKE) -C $(HTSSRC)
-
-
-
-#modied from htslib makefile
-FLAGS = -O3
-CPPFLAGS := $(filter-out -DNDEBUG,$(CPPFLAGS))
-FLAGS2 = $(CPPFLAGS) $(FLAGS) $(LDFLAGS)
-
-CFLAGS := $(FLAGS2) $(CFLAGS)
-CXXFLAGS := $(FLAGS2) $(CXXFLAGS)
+ifdef EXTRA_LIB
+LDFLAGS  := $(addprefix -L,$(EXTRA_LIB)) $(LDFLAGS)
+endif
 
-OPT=-Wl,-O2
-# filter -O2 optimization flag CXXFLAGS
-CXXFLAGS := $(filter-out -O2,$(CXXFLAGS))
-CXXFLAGS := $(subst $(OPT),,$(CXXFLAGS))
+ifdef EXTRA_LIBS
+LIBS += $(EXTRA_LIBS)
+endif
 
-CSRC = $(wildcard *.c)
-CXXSRC = $(wildcard *.cpp)
-OBJ = $(CSRC:.c=.o) $(CXXSRC:.cpp=.o)
+# --- Crypto library detektion ---
+HAVE_CRYPTO := $(shell echo 'int main(){}'|$(CXX) -x c++ - -lcrypto -o /dev/null 2>/dev/null && echo 0 || echo 1)
+LIBS := -lz -lm -lbz2 -llzma -lpthread -lcurl
 
-prefix      = /usr/local
-exec_prefix = $(prefix)
-bindir      = $(exec_prefix)/bin
+# --- Mål og bygning ---
+PROGRAM = metaDMG-cpp
+all: version.h $(PROGRAM) misc
 
-INSTALL = install
-INSTALL_DIR = $(INSTALL) -dm0755
-INSTALL_PROGRAM = $(INSTALL) -m0755
+ifeq ($(HAVE_CRYPTO),0)
+  $(info Crypto library is available to link; adding -lcrypto)
+  LIBS += -lcrypto
+else
+  $(info Crypto library is not available to link; skipping -lcrypto)
+endif
 
-PROGRAMS = metaDMG-cpp
+# --- HTSLIB håndtering ---
+# HTSSRC is an absolute path (e.g., $(CURDIR)/htslib or user-specified)
 
-all: $(PROGRAMS) misc
+ifndef HTSSRC
+  $(info HTSSRC not defined; cloning htslib from GitHub)
+  HTSSRC := $(CURDIR)/htslib
+  ABSPATH=$(HTSSRC) #donkykong
+endif
 
-PACKAGE_VERSION  = 0.4
+ABSPATH=$(HTSSRC) #donkykong
+ifeq ($(HTSSRC),systemwide)
+  $(info HTSSRC set to systemwide; using systemwide installation)
+  LIBS += -lhts
+  LIBHTS :=
+else
+  # Use HTSSRC directly for include path
+  CPPFLAGS += -I$(HTSSRC)
+  LIBHTS := $(HTSSRC)/libhts.a
+  LIBS := $(LIBHTS) $(LIBS)
+  $(PROGRAM): $(LIBHTS)
+
+  ifneq ($(filter /%,$(HTSSRC)),$(HTSSRC))
+    ABSPATH=../$(HTSSRC)
+  endif
+endif
 
-ifneq "$(wildcard .git)" ""
-PACKAGE_VERSION := $(shell git describe --always --dirty)
-version.h: $(if $(wildcard version.h),$(if $(findstring "$(PACKAGE_VERSION)",$(shell cat version.h)),,force))
+# Ensure htslib is cloned and built only if libhts.a is missing
+$(LIBHTS): .clone_htslib
+
+.clone_htslib:
+	@if [ ! -d "$(HTSSRC)" ]; then \
+		echo "Cloning htslib into $(HTSSRC) with submodules..."; \
+		git clone --recursive https://github.com/samtools/htslib.git $(HTSSRC) || { echo "Clone failed"; exit 1; }; \
+	else \
+		echo "$(HTSSRC) already exists, skipping clone."; \
+	fi
+	@if [ ! -f "$(LIBHTS)" ]; then \
+		$(MAKE) -C $(HTSSRC) libhts.a || { echo "Failed to build libhts.a"; exit 1; }; \
+	else \
+		echo "$(LIBHTS) already exists, skipping build."; \
+	fi
+	$(info CPPFLAGS=$(CPPFLAGS) HTSSRC=$(HTSSRC) (absolute path))
+
+# --- Versionsnummer og version.h ---
+PACKAGE_VERSION := 0.4
+ifneq ("$(wildcard .git)","")
+  PACKAGE_VERSION := $(shell git describe --always --dirty)
 endif
 
 version.h:
-	echo '#define METADAMAGE_VERSION "$(PACKAGE_VERSION)"' > $@
+	@echo '#define METADAMAGE_VERSION "$(PACKAGE_VERSION)"' > version.h.tmp
+	@cmp -s version.h.tmp version.h || mv version.h.tmp version.h
+	@rm -f version.h.tmp
 
-.PHONY: all clean install install-all install-misc misc test
+$(PROGRAM): version.h $(OBJ) $(LIBHTS)
+	$(CXX) $(CXXFLAGS) $(CPPFLAGS) -o $@ $(OBJ) $(LIBS) $(LDFLAGS)
 
-misc:
-	$(MAKE) -C misc HTSSRC="$(realpath $(HTSSRC))"
+.PHONY: misc
+misc: $(LIBHTS) $(OBJ)
+	$(MAKE) -C misc HTSSRC=$(ABSPATH)
 
+# --- Automatisk afhængighedshåndtering ---
 -include $(OBJ:.o=.d)
 
-%.o: %.c
-	$(CC) -c  $(CFLAGS) $*.c
-	$(CC) -MM $(CFLAGS) $*.c >$*.d
-
-%.o: %.cpp
-	$(CXX) -c  $(CXXFLAGS) $*.cpp
-	$(CXX) -MM $(CXXFLAGS) $*.cpp >$*.d
+%.o: %.c $(LIBHTS)
+	$(CC) -c $(CPPFLAGS) $(CFLAGS) $< -o $@
+	$(CC) -MM $(CPPFLAGS) $(CFLAGS) $< > $*.d
 
+%.o: %.cpp $(LIBHTS)
+	$(CXX) -c $(CPPFLAGS) $(CXXFLAGS) $< -o $@
+	$(CXX) -MM $(CPPFLAGS) $(CXXFLAGS) $< > $*.d
 
-metaDMG-cpp: version.h $(OBJ)
-	$(CXX) $(FLAGS) -o metaDMG-cpp *.o $(LIBS)
+# --- Rens og tests ---
+.PHONY: clean test testclean force
 
+clean: 
+	rm -f *.o *.d $(PROGRAM) version.h *~
+	rm -rf $(HTSSRC)
+	$(MAKE) -C misc clean
 
 testclean:
-	rm -f test/acc2taxid.map.gzf570b1db7c.dedup.filtered.rname.bam.bin
-	rm -f test/data/f570b1db7c.dedup.filtered.rname.bam
+	rm -f test/*.bam.bin test/data/*.bam
 	rm -rf test/output test/logfile version.h
 
-clean: testclean
-	rm  -f *.o *.d $(PROGRAMS) version.h *~
-	$(MAKE) -C misc clean
-
 test:
-	echo "Running unit tests for metaDMG"
-	cd test;./testAll.sh
+	@echo "Running unit tests for metaDMG"
+	cd test ; ./testAll.sh
 
 force:
-
-install: all
-	$(INSTALL_DIR) $(DESTDIR)$(bindir)
-	$(INSTALL_PROGRAM) $(PROGRAMS) $(DESTDIR)$(bindir)
-	$(MAKE) -C misc HTSSRC="$(realpath $(HTSSRC))" install
-
-install-misc: misc
-	$(MAKE) -C misc HTSSRC="$(realpath $(HTSSRC))" install-misc
-
-install-all: install install-misc

main_aggregate_stat.cpp

@@ -66,11 +66,11 @@ void to_root(int from,int to,std::map<int,mydata2> &stats,int2int &parent,int nr
     mydata2 &md2 = stats.find(to)->second;
     md2.data[0] = ((double) md1.data[0]*md1.nreads+md2.data[0]*md2.nreads)/((double) md1.nreads+md2.nreads); //weighted mean of read length
     md2.data[2] = ((double) md1.data[2]*md1.nreads+md2.data[2]*md2.nreads)/((double) md1.nreads+md2.nreads); //weighted mean of gc
-    
-    double variance1;
-    int nreads1;
-    double variance2;
-    int nreads2;
+
+#if 0
+    double variance1,variance2;
+    int nreads1, nreads2;
+#endif
     if(((double) md1.nreads+md2.nreads)>2){//pooled variance of length and GC
       /*variance1 = ((double) md1.nreads-1)*md1.data[1];
       nreads1 += md1.nreads;
@@ -130,7 +130,7 @@ std::map<int,char *> read_dfit(char *fname){
 }
 
 
-void aggr_stat2000(std::map<int, mydata2> &stats,int2int &parent){
+void aggr_stat3000(std::map<int, mydata2> &stats,int2int &parent){
   std::map<int,int> dasmap;
   for(std::map<int,mydata2>::iterator it = stats.begin();it!=stats.end();it++)
     dasmap[it->first] = it->second.nreads;
@@ -217,14 +217,6 @@ int main_aggregate(int argc, char **argv) {
     float postsize = retmap.size();
     fprintf(stderr, "\t-> pre: %f post:%f grownbyfactor: %f\n", presize, postsize, postsize / presize);
 
-    for (std::map<int, mydataD>::iterator it = retmap.begin(); it != retmap.end(); it++) {
-        int taxid = it->first;
-        mydataD md = it->second;
-        if (it->second.nal == 0)
-	  continue;
-	//	fprintf(stderr,"retmap taxid:%d nreads: %d\n",it->first,it->second.nreads);
-    }
-
     std::map<int, mydata2> stats;
     if (infile_lcastat){
       stats = load_lcastat(infile_lcastat,1);
@@ -248,7 +240,7 @@ int main_aggregate(int argc, char **argv) {
     if(child.size()>0)//this will not work if we have data at internal nodes 
       getval_stats(stats, child, 1);  // this will do everything
 #endif
-    aggr_stat2000(stats,parent);
+    aggr_stat3000(stats,parent);
     //    exit(0);
     if(stats.size()>0){
       postsize = stats.size();

merge_bdamage.cpp

@@ -89,7 +89,7 @@ void merge_bdamage(const std::vector<std::string> &bdamage_files, const char* ou
 void merge_rlens(const std::vector<std::string> &rlens_files, const char* outname) {
 
   // results 
-  std::map<int, std::array<size_t, 200>> rlens_merged;
+  std::map<int, std::array<size_t, 500>> rlens_merged;
 
   // each file 
   for (const std::string &fname : rlens_files) {
@@ -109,23 +109,23 @@ void merge_rlens(const std::vector<std::string> &rlens_files, const char* outnam
 
     while ((ret = bgzf_getline(fp, '\n', &line)) >= 0) {
       int taxid;
-      size_t tmp[200];
+      size_t tmp[500];
       char *ptr = line.s;
       char *endptr;
 
       // taxid is first column
       taxid = strtol(ptr, &endptr, 10); 
       ptr = endptr;
 
-      // each new cell (200)
-      for (int j = 0; j < 200; j++) {
+      // each new cell (500)
+      for (int j = 0; j < 500; j++) {
         tmp[j] = strtoull(ptr, &endptr, 10);
         ptr = endptr;
       }
 
       // sum and merge
       auto &entry = rlens_merged[taxid];
-      for (int j = 0; j < 200; j++) {
+      for (int j = 0; j < 500; j++) {
         entry[j] += tmp[j];
       }
     }
@@ -143,15 +143,15 @@ void merge_rlens(const std::vector<std::string> &rlens_files, const char* outnam
 
   kstring_t kstr = {0, 0, NULL};
   ksprintf(&kstr, "id");
-  for (int i = 0; i < 200; i++)
+  for (int i = 0; i < 500; i++)
     ksprintf(&kstr, "\trlen%d", i);
   ksprintf(&kstr, "\n");
 
   for (const auto &it : rlens_merged) {
     ksprintf(&kstr, "%d", it.first);
-    for (int i = 0; i < 199; i++)
+    for (int i = 0; i < 499; i++)
       ksprintf(&kstr, "\t%lu", it.second[i]);
-    ksprintf(&kstr, "\t%lu\n", it.second[199]);
+    ksprintf(&kstr, "\t%lu\n", it.second[499]);
 
     if (kstr.l > 1000000) {
       assert(bgzf_write(fp, kstr.s, kstr.l) == (ssize_t)kstr.l);