@@ -288,7 +288,7 @@ void genotype_hla(std::string ref_path,
288288
289289 {
290290 // Fix sequences names
291- auto const is_pass_alt = var.generate_infos ();
291+ auto is_pass_alt = var.generate_infos ();
292292 assert (var.seqs .size () == hla.sample_names .size ());
293293 assert ((is_pass_alt.size () + 1u ) == hla.sample_names .size ());
294294 std::vector<std::vector<char >> new_seqs;
@@ -306,6 +306,20 @@ void genotype_hla(std::string ref_path,
306306 new_seqs.push_back (seq);
307307 }
308308
309+ if (new_seqs.size () == 1 && hla.sample_names .size () >= 2 )
310+ {
311+ // special case where only reference allele is called, in this case also add some other sequence
312+ is_pass_alt[0 ] = true ;
313+ std::string const & hla_allele = hla.sample_names [1 ];
314+ std::vector<char > seq;
315+ seq.push_back (' <' );
316+ std::copy (hla_allele.begin (), hla_allele.end (), std::back_inserter (seq));
317+ seq.push_back (' >' );
318+ new_seqs.push_back (seq);
319+ }
320+
321+ assert (new_seqs.size () >= 2 );
322+
309323 for (auto const & new_seq : new_seqs)
310324 {
311325 print_log (log_severity::info, " Called contig sequence: " , std::string (new_seq.begin (), new_seq.end ()));
@@ -351,14 +365,15 @@ void genotype_hla(std::string ref_path,
351365
352366 var.seqs = std::move (new_seqs);
353367 var.stats = VarStats ();
368+ var.infos .clear ();
354369 var.generate_infos ();
355370 }
356371
357372 hla_vcf.open_for_writing (Options::const_instance ()->threads );
358373 hla_vcf.write_header ();
359374 assert (hla_vcf.variants .size () > 0 );
360375
361- if (hla_vcf. variants [ 0 ] .seqs .size () <= MAX_ALLELES)
376+ if (var .seqs .size () <= MAX_ALLELES)
362377 {
363378 hla_vcf.write_record (hla_vcf.variants [0 ], " .all" , false , false );
364379 }
0 commit comments