11package org .labkey .sequenceanalysis .run .alignment ;
22
3+ import htsjdk .samtools .util .CloseableIterator ;
4+ import htsjdk .variant .variantcontext .VariantContext ;
5+ import htsjdk .variant .vcf .VCFFileReader ;
36import org .json .JSONObject ;
47import org .labkey .api .module .ModuleLoader ;
58import org .labkey .api .pipeline .PipelineJob ;
1518import org .labkey .api .sequenceanalysis .pipeline .ToolParameterDescriptor ;
1619import org .labkey .api .sequenceanalysis .run .AbstractCommandWrapper ;
1720import org .labkey .api .sequenceanalysis .run .SimpleScriptWrapper ;
21+ import org .labkey .api .writer .PrintWriters ;
1822import org .labkey .sequenceanalysis .SequenceAnalysisModule ;
1923import org .labkey .sequenceanalysis .util .SequenceUtil ;
2024
2125import java .io .File ;
2226import java .io .IOException ;
27+ import java .io .PrintWriter ;
2328import java .util .ArrayList ;
2429import java .util .Arrays ;
2530import java .util .List ;
@@ -125,6 +130,12 @@ else if (!svVcf.exists())
125130 args .add (threads .toString ());
126131 }
127132
133+ args .add ("--region_file" );
134+ File regionFile = new File (ctx .getWorkingDirectory (), "regions.txt" );
135+ args .add (regionFile .getPath ());
136+
137+ generateRegionFile (svVcf , regionFile );
138+ ctx .getFileManager ().addIntermediateFile (regionFile );
128139 wrapper .execute (args );
129140
130141 File genotypes = new File (ctx .getWorkingDirectory (), "sv_results/" + SequenceAnalysisService .get ().getUnzippedBaseName (so .getName ()) + ".vcf.gz" );
@@ -146,4 +157,47 @@ else if (!svVcf.exists())
146157 }
147158 }
148159 }
160+
161+ private static void generateRegionFile (File svVcf , File regionFile ) throws PipelineJobException
162+ {
163+ try (PrintWriter writer = PrintWriters .getPrintWriter (regionFile ))
164+ {
165+ try (VCFFileReader reader = new VCFFileReader (svVcf , true ); CloseableIterator <VariantContext > iterator = reader .iterator ())
166+ {
167+ String chr = null ;
168+ int start = 1 ;
169+ int end = -1 ;
170+ while (iterator .hasNext ())
171+ {
172+ VariantContext vc = iterator .next ();
173+ if (chr == null )
174+ {
175+ chr = vc .getContig ();
176+ }
177+
178+ if (!vc .getContig ().equals (chr ))
179+ {
180+ writer .println (chr + ":" + start + "-" + end );
181+
182+ // Reset
183+ chr = vc .getContig ();
184+ start = 1 ;
185+ end = -1 ;
186+ }
187+
188+ start = Math .min (start , vc .getStart ());
189+ end = Math .max (end , vc .getEnd ());
190+ }
191+
192+ if (chr != null )
193+ {
194+ writer .println (chr + ":" + start + "-" + end );
195+ }
196+ }
197+ }
198+ catch (IOException e )
199+ {
200+ throw new PipelineJobException (e );
201+ }
202+ }
149203}
0 commit comments