Option to output full list of unplaced tags.
[lightning.git] / import.go
index 14b3667ed7fec34df0c75a4ef5504b4dc356d838..57596fbb07fb96b951ead13f41697e5f565e13d2 100644 (file)
--- a/import.go
+++ b/import.go
@@ -8,7 +8,6 @@ import (
        "flag"
        "fmt"
        "io"
-       "log"
        "net/http"
        _ "net/http/pprof"
        "os"
@@ -23,6 +22,7 @@ import (
        "time"
 
        "git.arvados.org/arvados.git/sdk/go/arvados"
+       log "github.com/sirupsen/logrus"
 )
 
 type importer struct {
@@ -31,6 +31,9 @@ type importer struct {
        outputFile     string
        projectUUID    string
        runLocal       bool
+       skipOOO        bool
+       outputTiles    bool
+       includeNoCalls bool
        encoder        *gob.Encoder
 }
 
@@ -45,10 +48,15 @@ func (cmd *importer) RunCommand(prog string, args []string, stdin io.Reader, std
        flags.SetOutput(stderr)
        flags.StringVar(&cmd.tagLibraryFile, "tag-library", "", "tag library fasta `file`")
        flags.StringVar(&cmd.refFile, "ref", "", "reference fasta `file`")
-       flags.StringVar(&cmd.outputFile, "o", "", "output `file`")
+       flags.StringVar(&cmd.outputFile, "o", "-", "output `file`")
        flags.StringVar(&cmd.projectUUID, "project", "", "project `UUID` for output data")
        flags.BoolVar(&cmd.runLocal, "local", false, "run on local host (default: run in an arvados container)")
+       flags.BoolVar(&cmd.skipOOO, "skip-ooo", false, "skip out-of-order tags")
+       flags.BoolVar(&cmd.outputTiles, "output-tiles", false, "include tile variant sequences in output file")
+       flags.BoolVar(&cmd.includeNoCalls, "include-no-calls", false, "treat tiles with no-calls as regular tiles")
+       priority := flags.Int("priority", 500, "container request priority")
        pprof := flags.String("pprof", "", "serve Go profile data at http://`[addr]:port`")
+       loglevel := flags.String("loglevel", "info", "logging threshold (trace, debug, info, warn, error, fatal, or panic)")
        err = flags.Parse(args)
        if err == flag.ErrHelp {
                err = nil
@@ -69,13 +77,20 @@ func (cmd *importer) RunCommand(prog string, args []string, stdin io.Reader, std
                }()
        }
 
+       lvl, err := log.ParseLevel(*loglevel)
+       if err != nil {
+               return 2
+       }
+       log.SetLevel(lvl)
+
        if !cmd.runLocal {
                runner := arvadosContainerRunner{
                        Name:        "lightning import",
                        Client:      arvados.NewClientFromEnv(),
                        ProjectUUID: cmd.projectUUID,
-                       RAM:         30000000000,
+                       RAM:         60000000000,
                        VCPUs:       16,
+                       Priority:    *priority,
                }
                err = runner.TranslatePaths(&cmd.tagLibraryFile, &cmd.refFile, &cmd.outputFile)
                if err != nil {
@@ -88,7 +103,7 @@ func (cmd *importer) RunCommand(prog string, args []string, stdin io.Reader, std
                                return 1
                        }
                }
-               if cmd.outputFile == "" {
+               if cmd.outputFile == "-" {
                        cmd.outputFile = "/mnt/output/library.gob"
                } else {
                        // Not yet implemented, but this should write
@@ -97,11 +112,22 @@ func (cmd *importer) RunCommand(prog string, args []string, stdin io.Reader, std
                        err = errors.New("cannot specify output file in container mode: not implemented")
                        return 1
                }
-               runner.Args = append([]string{"import", "-local=true", "-tag-library", cmd.tagLibraryFile, "-ref", cmd.refFile, "-o", cmd.outputFile}, inputs...)
-               err = runner.Run()
+               runner.Args = append([]string{"import",
+                       "-local=true",
+                       "-loglevel=" + *loglevel,
+                       fmt.Sprintf("-skip-ooo=%v", cmd.skipOOO),
+                       fmt.Sprintf("-output-tiles=%v", cmd.outputTiles),
+                       fmt.Sprintf("-include-no-calls=%v", cmd.includeNoCalls),
+                       "-tag-library", cmd.tagLibraryFile,
+                       "-ref", cmd.refFile,
+                       "-o", cmd.outputFile,
+               }, inputs...)
+               var output string
+               output, err = runner.Run()
                if err != nil {
                        return 1
                }
+               fmt.Fprintln(stdout, output+"/library.gob")
                return 0
        }
 
@@ -110,18 +136,13 @@ func (cmd *importer) RunCommand(prog string, args []string, stdin io.Reader, std
                return 1
        }
 
-       tilelib, err := cmd.loadTileLibrary()
+       taglib, err := cmd.loadTagLibrary()
        if err != nil {
                return 1
        }
-       go func() {
-               for range time.Tick(10 * time.Second) {
-                       log.Printf("tilelib.Len() == %d", tilelib.Len())
-               }
-       }()
 
        var output io.WriteCloser
-       if cmd.outputFile == "" {
+       if cmd.outputFile == "-" {
                output = nopCloser{stdout}
        } else {
                output, err = os.OpenFile(cmd.outputFile, os.O_CREATE|os.O_WRONLY, 0777)
@@ -133,6 +154,17 @@ func (cmd *importer) RunCommand(prog string, args []string, stdin io.Reader, std
        bufw := bufio.NewWriter(output)
        cmd.encoder = gob.NewEncoder(bufw)
 
+       tilelib := &tileLibrary{taglib: taglib, includeNoCalls: cmd.includeNoCalls, skipOOO: cmd.skipOOO}
+       if cmd.outputTiles {
+               cmd.encoder.Encode(LibraryEntry{TagSet: taglib.Tags()})
+               tilelib.encoder = cmd.encoder
+       }
+       go func() {
+               for range time.Tick(10 * time.Minute) {
+                       log.Printf("tilelib.Len() == %d", tilelib.Len())
+               }
+       }()
+
        err = cmd.tileInputs(tilelib, infiles)
        if err != nil {
                return 1
@@ -165,7 +197,7 @@ func (cmd *importer) tileFasta(tilelib *tileLibrary, infile string) (tileSeq, er
        return tilelib.TileFasta(infile, input)
 }
 
-func (cmd *importer) loadTileLibrary() (*tileLibrary, error) {
+func (cmd *importer) loadTagLibrary() (*tagLibrary, error) {
        log.Printf("tag library %s load starting", cmd.tagLibraryFile)
        f, err := os.Open(cmd.tagLibraryFile)
        if err != nil {
@@ -189,15 +221,22 @@ func (cmd *importer) loadTileLibrary() (*tileLibrary, error) {
                return nil, fmt.Errorf("cannot tile: tag library is empty")
        }
        log.Printf("tag library %s load done", cmd.tagLibraryFile)
-       return &tileLibrary{taglib: &taglib}, nil
+       return &taglib, nil
 }
 
+var (
+       vcfFilenameRe    = regexp.MustCompile(`\.vcf(\.gz)?$`)
+       fasta1FilenameRe = regexp.MustCompile(`\.1\.fa(sta)?(\.gz)?$`)
+       fasta2FilenameRe = regexp.MustCompile(`\.2\.fa(sta)?(\.gz)?$`)
+       fastaFilenameRe  = regexp.MustCompile(`\.fa(sta)?(\.gz)?$`)
+)
+
 func listInputFiles(paths []string) (files []string, err error) {
        for _, path := range paths {
                if fi, err := os.Stat(path); err != nil {
                        return nil, fmt.Errorf("%s: stat failed: %s", path, err)
                } else if !fi.IsDir() {
-                       if !strings.HasSuffix(path, ".2.fasta") || strings.HasSuffix(path, ".2.fasta.gz") {
+                       if !fasta2FilenameRe.MatchString(path) {
                                files = append(files, path)
                        }
                        continue
@@ -213,23 +252,27 @@ func listInputFiles(paths []string) (files []string, err error) {
                }
                sort.Strings(names)
                for _, name := range names {
-                       if strings.HasSuffix(name, ".vcf") || strings.HasSuffix(name, ".vcf.gz") {
+                       if vcfFilenameRe.MatchString(name) {
                                files = append(files, filepath.Join(path, name))
-                       } else if strings.HasSuffix(name, ".1.fasta") || strings.HasSuffix(name, ".1.fasta.gz") {
+                       } else if fastaFilenameRe.MatchString(name) && !fasta2FilenameRe.MatchString(name) {
                                files = append(files, filepath.Join(path, name))
                        }
                }
                d.Close()
        }
        for _, file := range files {
-               if strings.HasSuffix(file, ".1.fasta") || strings.HasSuffix(file, ".1.fasta.gz") {
-                       continue
-               } else if _, err := os.Stat(file + ".csi"); err == nil {
-                       continue
-               } else if _, err = os.Stat(file + ".tbi"); err == nil {
+               if fastaFilenameRe.MatchString(file) {
                        continue
+               } else if vcfFilenameRe.MatchString(file) {
+                       if _, err := os.Stat(file + ".csi"); err == nil {
+                               continue
+                       } else if _, err = os.Stat(file + ".tbi"); err == nil {
+                               continue
+                       } else {
+                               return nil, fmt.Errorf("%s: cannot read without .tbi or .csi index file", file)
+                       }
                } else {
-                       return nil, fmt.Errorf("%s: cannot read without .tbi or .csi index file", file)
+                       return nil, fmt.Errorf("don't know how to handle filename %s", file)
                }
        }
        return
@@ -245,25 +288,42 @@ func (cmd *importer) tileInputs(tilelib *tileLibrary, infiles []string) error {
                var phases sync.WaitGroup
                phases.Add(2)
                variants := make([][]tileVariantID, 2)
-               if strings.HasSuffix(infile, ".1.fasta") || strings.HasSuffix(infile, ".1.fasta.gz") {
+               if fasta1FilenameRe.MatchString(infile) {
                        todo <- func() error {
                                defer phases.Done()
                                log.Printf("%s starting", infile)
                                defer log.Printf("%s done", infile)
                                tseqs, err := cmd.tileFasta(tilelib, infile)
-                               variants[0] = tseqs.Variants()
+                               var kept, dropped int
+                               variants[0], kept, dropped = tseqs.Variants()
+                               log.Printf("%s found %d unique tags plus %d repeats", infile, kept, dropped)
                                return err
                        }
-                       infile2 := regexp.MustCompile(`\.1\.fasta(\.gz)?$`).ReplaceAllString(infile, `.2.fasta$1`)
+                       infile2 := fasta1FilenameRe.ReplaceAllString(infile, `.2.fa$1$2`)
                        todo <- func() error {
                                defer phases.Done()
                                log.Printf("%s starting", infile2)
                                defer log.Printf("%s done", infile2)
                                tseqs, err := cmd.tileFasta(tilelib, infile2)
-                               variants[1] = tseqs.Variants()
+                               var kept, dropped int
+                               variants[1], kept, dropped = tseqs.Variants()
+                               log.Printf("%s found %d unique tags plus %d repeats", infile, kept, dropped)
                                return err
                        }
-               } else {
+               } else if fastaFilenameRe.MatchString(infile) {
+                       todo <- func() error {
+                               defer phases.Done()
+                               defer phases.Done()
+                               log.Printf("%s starting", infile)
+                               defer log.Printf("%s done", infile)
+                               tseqs, err := cmd.tileFasta(tilelib, infile)
+                               var kept, dropped int
+                               variants[0], kept, dropped = tseqs.Variants()
+                               variants[1] = variants[0]
+                               log.Printf("%s found %d unique tags plus %d repeats", infile, kept, dropped)
+                               return err
+                       }
+               } else if vcfFilenameRe.MatchString(infile) {
                        for phase := 0; phase < 2; phase++ {
                                phase := phase
                                todo <- func() error {
@@ -271,10 +331,14 @@ func (cmd *importer) tileInputs(tilelib *tileLibrary, infiles []string) error {
                                        log.Printf("%s phase %d starting", infile, phase+1)
                                        defer log.Printf("%s phase %d done", infile, phase+1)
                                        tseqs, err := cmd.tileGVCF(tilelib, infile, phase)
-                                       variants[phase] = tseqs.Variants()
+                                       var kept, dropped int
+                                       variants[phase], kept, dropped = tseqs.Variants()
+                                       log.Printf("%s phase %d found %d unique tags plus %d repeats", infile, phase+1, kept, dropped)
                                        return err
                                }
                        }
+               } else {
+                       panic(fmt.Sprintf("bug: unhandled filename %q", infile))
                }
                encodeJobs.Add(1)
                go func() {
@@ -289,8 +353,11 @@ func (cmd *importer) tileInputs(tilelib *tileLibrary, infiles []string) error {
                        }
                        flat := make([]tileVariantID, ntags*2)
                        for i := 0; i < ntags; i++ {
-                               flat[i*2] = variants[0][i]
-                               flat[i*2+1] = variants[1][i]
+                               for hap := 0; hap < 2; hap++ {
+                                       if i < len(variants[hap]) {
+                                               flat[i*2+hap] = variants[hap][i]
+                                       }
+                               }
                        }
                        err := cmd.encoder.Encode(LibraryEntry{
                                CompactGenomes: []CompactGenome{{Name: infile, Variants: flat}},
@@ -305,9 +372,10 @@ func (cmd *importer) tileInputs(tilelib *tileLibrary, infiles []string) error {
        }
        go close(todo)
        var tileJobs sync.WaitGroup
-       running := int64(runtime.NumCPU())
-       for i := 0; i < runtime.NumCPU(); i++ {
+       var running int64
+       for i := 0; i < runtime.NumCPU()*9/8+1; i++ {
                tileJobs.Add(1)
+               atomic.AddInt64(&running, 1)
                go func() {
                        defer tileJobs.Done()
                        defer atomic.AddInt64(&running, -1)