diff CSP2/CSP2_env/env-d9b9114564458d9d-741b3de822f2aaca6c6caa4325c4afce/opt/mummer-3.23/nucmer2xfig @ 69:33d812a61356

planemo upload commit 2e9511a184a1ca667c7be0c6321a36dc4e3d116d
author jpayne
date Tue, 18 Mar 2025 17:55:14 -0400
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/CSP2/CSP2_env/env-d9b9114564458d9d-741b3de822f2aaca6c6caa4325c4afce/opt/mummer-3.23/nucmer2xfig	Tue Mar 18 17:55:14 2025 -0400
@@ -0,0 +1,139 @@
+#!/usr/bin/env perl
+# (c) Steven Salzberg 2001
+# Make an xfig plot for a comparison of a reference chromosome (or single
+# molecule) versus a multifasta file of contigs from another genome.  
+# The input file here is a NUCmer coords file such as:
+# 2551  2577 | 240 266 | 27  27 | 96.30  | 20302755  1424  | 0.00  1.90 | 2R 1972084
+# generated by running 'show-coords -c -l <delta file>'
+# For the above example, D. melanogaster chr 2R is the reference and the query
+# is an assembly with 1000s of contigs from D. pseudoobscura
+# The file needs to be sorted by the smaller contig ids, via:
+# tail +6 Dmel2R-vs-Dpseudo.coords | sort -k 19n -k 4n  > Dmel2R-vs-Dpseudo-resort.coords 
+#  and remember to get rid of top 5 (header) lines.
+# Usage: plot-drosoph-align-xfig.perl Dmel2R-vs-Dpseudo-resort.coords 
+unless (open(coordsfile,$ARGV[0]))  {
+    die ("can't open file $ARGV[0].\n");
+}
+$Xscale = 0.005;
+$Yscale = 20;
+$chrcolor = 4;   # 4 is red, 1 is blue, 2 is green
+$green = 2;
+$contigcolor = 1; # contigs are blue
+
+# print header info
+print "#FIG 3.2\nLandscape\nCenter\nInches\nLetter  \n100.00\nSingle\n-2\n1200 2\n";
+$first_time = 1;
+while (<coordsfile>) {
+    ($s1,$e1,$x1,$s2,$e2,$x2,$l1,$l2,$x3,$percent_id,$x4,$Rlen,$Qlen,$x5,$Rcov,$Qcov,$x6,$Rid,$Qid) = split(" ");
+    if ($prevQid eq $Qid) {  # query contig is same as prev line
+	$dist = abs($s1 - $prev_s1);
+	if ( $dist >  2 * $Qlen)  { # if this match is too far away
+	    # print the contig here if the matching bit is > 1000
+	    if ($right_end{$Qid} - $left_end{$Qid} > 1000) {
+		# print it at y=50 rather than 100 because the scale is 0-50,
+		# where 0=50% identical and 50 is 100% identical
+		print_xfig_line($left_end{$Qid},50,$right_end{$Qid},50,$contigcolor);
+		print_label($left_end{$Qid},50,$right_end{$Qid},$Qid,5);
+	    }
+	    $left_end{$Qid} = $s1;    # then re-set the start and end of the contig
+	    $right_end{$Qid} = $e1;   # and we'll print it again later
+	}
+	else {   # extend the boundaries of the match
+	    if ($s1 < $left_end{$Qid}) { $left_end{$Qid} = $s1; }
+	    if ($e1 > $right_end{$Qid}) { $right_end{$Qid} = $e1; }
+	}
+    }
+    else { # this is a different contig, first time seeing it
+	$left_end{$Qid} = $s1;
+	$right_end{$Qid} = $e1;
+	# print the previous contig as a line at y=100 for 100%
+	if ($first_time < 1) {
+	    print_xfig_line($left_end{$prevQid},50,$right_end{$prevQid},50,$contigcolor);
+	    print_label($left_end{$prevQid},50,$right_end{$prevQid},$prevQid,5);
+	}
+	else { $first_time = 0; }
+    }
+    $prevQid = $Qid;
+    $prev_s1 = $s1;
+    $prev_Qlen = $Qlen;
+    # next print the matching bit as a separate line, with a separate color,
+    # with its height determined by percent match
+    $Xleft = int($Xscale * $s1);
+    $Xright = int($Xscale * $e1);
+    if ($Xleft == $Xright) { $Xright += 1; }
+    if ($percent_id < 50) { $percent_id = 50; }
+    print_xfig_line($s1,$percent_id - 50,$e1,$percent_id - 50,$green);
+}
+# print very last contig
+$left_end{$Qid} = $s1;
+$right_end{$Qid} = $e1;
+print_xfig_line($s1,50,$e1,50,$contigcolor);
+print_label($s1,50,$e1,$Qid,5);
+close(coordsfile);
+
+# now draw the horizontal chr line for the reference
+$label_xpos = $Xscale * $Rlen;
+print "4 0 0 100 0 0 12 0.0000 4 135 405 ";
+printf("%.0f %.0f",$label_xpos, 0);
+print " ", $Rid, "\\001\n";
+print "2 1 0 2 $chrcolor 7 100 0 -1 0.000 0 0 -1 0 0 2\n";		
+printf("\t %.0f %.0f %.0f %.0f\n", 0, 0,
+       $Xscale * $Rlen, 0);
+# print some X-axis coordinates
+$pointsize = 5;
+for ($i = 250000; $i < $Rlen - 50000; $i+= 250000) { 
+    print "4 0 0 100 0 0 $pointsize 0.0000 4 135 405 ";
+    printf("%.0f %.0f",$i * $Xscale + 20, -20 + ($Yscale * -0.5));
+    print " $i", "\\001\n";
+    #print a vertical tic mark
+    print "2 1 0 1 0 7 50 0 -1 0.000 0 0 -1 0 0 2\n"; 
+    printf("%.0f %.0f %.0f -50\n",
+	   $i * $Xscale, 0, $i * $Xscale);
+}
+# print tic marks indicating % identity on the y-axis
+for ($percent_id = 50; $percent_id < 101; $percent_id += 10) {
+    print "4 0 0 100 0 0 $pointsize 0.0000 4 135 405 ";
+    printf("-150 %.0f", ($percent_id - 50) * $Yscale + 10); # shift down 10 pixels
+    print " $percent_id", "\\001\n";
+    # print the tic mark
+    print "2 1 0 1 0 7 50 0 -1 0.000 0 0 -1 0 0 2\n"; 
+    printf("-50 %.0f 0 %.0f\n",
+	   ($percent_id - 50) * $Yscale, ($percent_id - 50) * $Yscale);
+}
+
+# print a line in the appropriate color, scaled with Xscale,Yscale
+sub print_xfig_line {
+    my ($xleft,$yleft,$xright,$yright,$color) = @_;
+    # print it at the given coordinates, scaled
+    $xleft_scaled = int($Xscale * $xleft);
+    $xright_scaled = int($Xscale * $xright);
+    # Xfig has a bug: if we re-scale and the resulting X coordinates are equal,
+    # then it will print a fixed-size (large) rectangle, rather than a 1-pixel
+    # wide one.  So check fo this and correct.
+    if ($xleft_scaled == $xright_scaled) { $xright_scaled += 1; }
+    # set up and print line in xfig format
+    print "2 1 0 2 $color 7 100 0 -1 0.000 0 0 -1 0 0 2\n";		
+    printf("\t %.0f %.0f %.0f %.0f\n", 
+	   $xleft_scaled, $Yscale * $yleft, $xright_scaled, $Yscale * $yright);
+}
+
+# print a label for each contig using its ID
+sub print_label {
+    my ($xleft,$yleft,$xright,$id,$psize) = @_;
+    # print it at the left edge of the contig.  The angle
+    # of 4.7124 means text goes down vertically.  5th argument
+    # here is pointsize of the label.  
+    # bump the label down 20 pixels
+    $yposition = $yleft * $Yscale + 20;
+    # to keep the display clean, don't print the label unless the
+    # contig itself is wider than the width of the text in the label
+    $xleft_scaled = $xleft * $Xscale;
+    $xright_scaled = $xright * $Xscale;
+    # each point of type is 8 pixels
+    $textheight = $psize * 8;
+    if ($xright_scaled - $xleft_scaled > $textheight) {
+	print "4 0 0 50 0 0 $psize 4.7124 4 135 435 ";
+	printf("%.0f %.0f",$xleft_scaled, $yposition);
+	print " $id", "\\001\n";
+    }
+}