Merge pull request #51 from ComparativeGenomicsToolkit/bigdel

PAF deletion filter
ComparativeGenomicsToolkit · Dec 20, 2021 · 56ee085 · 56ee085
2 parents bfa3de8 + ce52e53
commit 56ee085
Show file tree

Hide file tree

Showing 3 changed files with 619 additions and 5 deletions.
diff --git a/Makefile b/Makefile
@@ -2,7 +2,7 @@
 rootPath = ./
 include ./include.mk
 
-all : hal2vg clip-vg halRemoveDupes halMergeChroms halUnclip
+all : hal2vg clip-vg halRemoveDupes halMergeChroms halUnclip filter-paf-deletions
 
 # Note: hdf5 from apt doesn't seem to work for static builds.  It should be installed
 # from source and configured with "--enable-static --disable-shared", then have its
@@ -38,12 +38,18 @@ ifeq ($(shell ldd halUnclip | grep "not a dynamic" | wc -l), $(shell ls halUncli
 else
 	$(error ldd found dnymaic linked dependency in halUnclip)
 endif
+ifeq ($(shell ldd filter-paf-deletions | grep "not a dynamic" | wc -l), $(shell ls filter-paf-deletions | wc -l))
+	$(info ldd verified that filter-paf-deletions static)
+else
+	$(error ldd found dnymaic linked dependency in filter-paf-deletions)
+endif
+
 
 cleanFast : 
-	rm -f hal2vg hal2vg.o clip-vg clip-vg.o halRemoveDupes halRemoveDupes.o halMergeChroms halMergeChroms.o halUnclip halUnclip.o
+	rm -f hal2vg hal2vg.o clip-vg clip-vg.o halRemoveDupes halRemoveDupes.o halMergeChroms halMergeChroms.o halUnclip halUnclip.o filter-paf-deletions filter-paf-deletions.o
 
 clean :
-	rm -f hal2vg hal2vg.o clip-vg clip-vg.o halRemoveDupes halRemoveDupes.o halMergeChroms halMergeChroms.o halUnclip halUnclip.o
+	rm -f hal2vg hal2vg.o clip-vg clip-vg.o halRemoveDupes halRemoveDupes.o halMergeChroms halMergeChroms.o halUnclip halUnclip.o filter-paf-deletions filter-paf-deletions.o
 	cd deps/sonLib && make clean
 	cd deps/pinchesAndCacti && make clean
 	cd deps/hal && make clean
@@ -91,6 +97,12 @@ halUnclip.o : halUnclip.cpp subpaths.h ${basicLibsDependencies}
 halUnclip : halUnclip.o ${basicLibsDependencies} 
 	${cpp} ${CXXFLAGS} -fopenmp -pthread halUnclip.o  ${basicLibs}  -o halUnclip
 
+filter-paf-deletions.o : filter-paf-deletions.cpp subpaths.h ${basicLibsDependencies}
+	${cpp} ${CXXFLAGS} -I . filter-paf-deletions.cpp -c
+
+filter-paf-deletions : filter-paf-deletions.o ${basicLibsDependencies} 
+	${cpp} ${CXXFLAGS} -fopenmp -pthread filter-paf-deletions.o  ${basicLibs}  -o filter-paf-deletions
+
 test :
 	make
 	cd tests && prove -v t
diff --git a/clip-vg.cpp b/clip-vg.cpp
@@ -768,11 +768,25 @@ void forwardize_paths(MutablePathMutableHandleGraph* graph, const string& ref_pr
                         if (graph->get_is_reverse(handle)) {
                             handle_t flipped_handle = graph->create_handle(graph->get_sequence(handle));
 			    graph->follow_edges(handle, true, [&](handle_t prev_handle) {
-				graph->create_edge(prev_handle, flipped_handle);
+                                    if (graph->get_id(prev_handle) != graph->get_id(handle)) {
+                                        graph->create_edge(prev_handle, flipped_handle);
+                                    }
 			      });
 			    graph->follow_edges(handle, false, [&](handle_t next_handle) {
-				graph->create_edge(flipped_handle, next_handle);
+                                    if (graph->get_id(handle) != graph->get_id(next_handle)) {
+                                        graph->create_edge(flipped_handle, next_handle);
+                                    }
 			      });
+                            // self-loop cases we punted on above:
+                            if (graph->has_edge(handle, handle)) {
+                                graph->create_edge(flipped_handle, flipped_handle);
+                            }
+                            if (graph->has_edge(handle, graph->flip(handle))) {
+                                graph->create_edge(flipped_handle, graph->flip(flipped_handle));                                
+                            }
+                            if (graph->has_edge(graph->flip(handle), handle)) {
+                                graph->create_edge(graph->flip(flipped_handle), flipped_handle);
+                            }
                             vector<step_handle_t> steps = graph->steps_of_handle(handle);
                             size_t ref_count = 0;
                             for (step_handle_t step : steps) {