afee0fcf2a7b7de2eb37e1420dc126607373453e hiram Wed May 27 17:31:45 2020 -0700 rework scripts to allow subset hubs to be created from a list of assemblies refs #23734 diff --git src/hg/makeDb/doc/asmHubs/trackData.pl src/hg/makeDb/doc/asmHubs/trackData.pl index 71fe9ff..ab7922d 100755 --- src/hg/makeDb/doc/asmHubs/trackData.pl +++ src/hg/makeDb/doc/asmHubs/trackData.pl @@ -1,53 +1,55 @@ #!/usr/bin/env perl use strict; use warnings; use File::stat; my $argc = scalar(@ARGV); -if ($argc < 2) { - printf STDERR "usage: trackData.pl Name asmHubName > trackData.html\n"; - printf STDERR "e.g.: trackData.pl Mammals mammals > trackData.html\n"; +if ($argc < 3) { + printf STDERR "usage: trackData.pl Name asmHubName [two column name list] > trackData.html\n"; + printf STDERR "e.g.: trackData.pl Mammals mammals mammals.asmId.commonName.tsv > trackData.html\n"; + printf STDERR "the name list is found in \$HOME/kent/src/hg/makeDb/doc/asmHubs/\n"; + printf STDERR "\nthe two columns are 1: asmId (accessionId_assemblyName)\n"; + printf STDERR "column 2: common name for species, columns separated by tab\n"; exit 255; } my $testOutput = 0; my $spliceOut = -1; if ($argc > 2) { for (my $i = 0; $i < $argc; ++$i) { if ($ARGV[$i] =~ /-test/) { $testOutput = 1; $spliceOut = $i; } } } if ($spliceOut != -1) { splice @ARGV, $spliceOut, 1; } my $Name = shift; my $asmHubName = shift; +my $commonNameOrder = shift; my $home = $ENV{'HOME'}; my $toolsDir = "$home/kent/src/hg/makeDb/doc/asmHubs"; -my $commonNameList = "$asmHubName.asmId.commonName.tsv"; -my $commonNameOrder = "$asmHubName.commonName.asmId.orderList.tsv"; -my @orderList; # asmId of the assemblies in order from the *.list files -# the order to read the different .list files: -my %betterName; # key is asmId, value is better common name than found in +my @orderList; # asmId of the assemblies in order from the commonNameOrder file +my %commonName; # key is asmId, value is a common name, perhaps more appropriate + # than found in assembly_report file # assembly_report my $vgpIndex = 0; $vgpIndex = 1 if ($Name =~ m/vgp/i); my $assemblyTotal = 0; # complete list of assemblies in this group my $asmCount = 0; # count of assemblies completed and in the table my $overallNucleotides = 0; my $overallSeqCount = 0; my $overallGapSize = 0; my $overallGapCount = 0; ############################################################################## # from Perl Cookbook Recipe 2.17, print out large numbers with comma delimiters: ############################################################################## sub commify($) { @@ -340,31 +342,31 @@ my $asmCounted = 0; if ($testOutput) { # add extra columns during 'test' output # 0 1 2 3 4 5 6 # 7 8 9 10 11 12 13 # 14 # my @trackList = qw(ncbiRefSeq xenoRefGene augustus ensGene gc5Base gap allGaps assembly rmsk simpleRepeat windowMasker gapOverlap tandemDups cpgIslandExtUnmasked cpgIslandExt); # 0 1 2 3 4 5 6 # 7 8 9 10 # my @trackList = qw(ncbiRefSeq xenoRefGene augustus ensGene gc5Base allGaps assembly rmsk simpleRepeat windowMasker cpgIslandExtUnmasked); splice @trackList, 11, 0, "cpgIslandExt"; splice @trackList, 10, 0, "tandemDups"; splice @trackList, 10, 0, "gapOverlap"; splice @trackList, 5, 0, "gap"; splice @trackList, 1, 0, "ncbiGene"; } - foreach my $asmId (reverse(@orderList)) { + foreach my $asmId (@orderList) { my $tracksCounted = 0; my ($gcPrefix, $asmAcc, $asmName) = split('_', $asmId, 3); my $accessionId = sprintf("%s_%s", $gcPrefix, $asmAcc); my $accessionDir = substr($asmId, 0 ,3); $accessionDir .= "/" . substr($asmId, 4 ,3); $accessionDir .= "/" . substr($asmId, 7 ,3); $accessionDir .= "/" . substr($asmId, 10 ,3); my $buildDir = "/hive/data/genomes/asmHubs/refseqBuild/$accessionDir/$asmId"; if ($gcPrefix eq "GCA") { $buildDir = "/hive/data/genomes/asmHubs/genbankBuild/$accessionDir/$asmId"; } my $trackDb="$buildDir/${asmId}.trackDb.txt"; next if (! -s "$trackDb"); # assembly build not complete my $asmReport="$buildDir/download/${asmId}_assembly_report.txt"; if (! -s "$asmReport") { @@ -405,31 +407,31 @@ $line =~ s/ //g;; $line =~ s/\s+$//g;; if ($line =~ m/Date:/) { if ($asmDate =~ m/notFound/) { ++$itemsFound; $asmDate = $line; $asmDate =~ s/.*:\s+//; } } elsif ($line =~ m/Organism name:/) { if ($sciName =~ m/notFound/) { ++$itemsFound; $commonName = $line; $sciName = $line; $commonName =~ s/.*\(//; $commonName =~ s/\)//; - $commonName = $betterName{$asmId} if (exists($betterName{$asmId})); + $commonName = $commonName{$asmId} if (exists($commonName{$asmId})); $sciName =~ s/.*:\s+//; $sciName =~ s/\s+\(.*//; } } } close (FH); my $hubUrl = "https://genome.ucsc.edu/h/$accessionId"; if ($testOutput) { $hubUrl = "https://genome-test.gi.ucsc.edu/h/$accessionId"; } printf "<tr><td align=right>%d</td>\n", ++$asmCount; printf "<td align=center><a href='%s' target=_blank>%s<br>%s</a></td>\n", $hubUrl, $commonName, $accessionId; foreach my $track (@trackList) { my $trackFile = "$buildDir/bbi/$asmId.$track"; my $trackFb = "$buildDir/trackData/$track/fb.$asmId.$track.txt"; @@ -493,27 +495,27 @@ } # foreach my $track (@trackList) printf "</tr>\n"; $asmCounted += 1; printf STDERR "# %03d\t%02d tracks\t%s\n", $asmCounted, $tracksCounted, $asmId; } } ############################################################################## ### main() ############################################################################## open (FH, "<$toolsDir/${commonNameOrder}") or die "can not read ${commonNameOrder}"; while (my $line = <FH>) { next if ($line =~ m/^#/); chomp $line; - my ($commonName, $asmId) = split('\t', $line); + my ($asmId, $commonName) = split('\t', $line); push @orderList, $asmId; - $betterName{$asmId} = $commonName; + $commonName{$asmId} = $commonName; ++$assemblyTotal; } close (FH); startHtml(); startTable(); tableContents(); endTable(); endHtml();