1 // Author: Derek Barnett
2 
3 #include <cstddef>
4 #include <string>
5 
6 #include <gtest/gtest.h>
7 
8 #include "PbbamTestData.h"
9 
10 #include <pbbam/EntireFileQuery.h>
11 #include <pbbam/PbiFilter.h>
12 #include <pbbam/virtual/VirtualPolymeraseCompositeReader.h>
13 #include <pbbam/virtual/VirtualPolymeraseReader.h>
14 #include <pbbam/virtual/ZmwReadStitcher.h>
15 
16 using namespace PacBio;
17 using namespace PacBio::BAM;
18 
19 namespace ZmwReadStitcherTests {
20 
Compare(const BamRecord & b1,const BamRecord & b2)21 static void Compare(const BamRecord& b1, const BamRecord& b2)
22 {
23     EXPECT_TRUE(b1.HasDeletionQV());
24     EXPECT_TRUE(b1.HasDeletionTag());
25     EXPECT_TRUE(b1.HasInsertionQV());
26     EXPECT_TRUE(b1.HasMergeQV());
27     EXPECT_TRUE(b1.HasSubstitutionQV());
28     EXPECT_TRUE(b1.HasSubstitutionTag());
29     EXPECT_TRUE(b1.HasLabelQV());
30     EXPECT_TRUE(b1.HasAltLabelQV());
31     EXPECT_TRUE(b1.HasAltLabelTag());
32     EXPECT_TRUE(b1.HasPkmean());
33     EXPECT_TRUE(b1.HasPkmid());
34     EXPECT_TRUE(b1.HasPulseCall());
35     EXPECT_TRUE(b1.HasIPD());
36     EXPECT_TRUE(b1.HasPulseWidth());
37     EXPECT_TRUE(b1.HasPrePulseFrames());
38     EXPECT_TRUE(b1.HasPulseCallWidth());
39     EXPECT_TRUE(b1.HasPulseMergeQV());
40 
41     EXPECT_TRUE(b2.HasDeletionQV());
42     EXPECT_TRUE(b2.HasDeletionTag());
43     EXPECT_TRUE(b2.HasInsertionQV());
44     EXPECT_TRUE(b2.HasMergeQV());
45     EXPECT_TRUE(b2.HasSubstitutionQV());
46     EXPECT_TRUE(b2.HasSubstitutionTag());
47     EXPECT_TRUE(b2.HasLabelQV());
48     EXPECT_TRUE(b2.HasAltLabelQV());
49     EXPECT_TRUE(b2.HasAltLabelTag());
50     EXPECT_TRUE(b2.HasPkmean());
51     EXPECT_TRUE(b2.HasPkmid());
52     EXPECT_TRUE(b2.HasPulseCall());
53     EXPECT_TRUE(b2.HasIPD());
54     EXPECT_TRUE(b2.HasPulseWidth());
55     EXPECT_TRUE(b2.HasPrePulseFrames());
56     EXPECT_TRUE(b2.HasPulseCallWidth());
57     EXPECT_TRUE(b2.HasPulseMergeQV());
58 
59     EXPECT_EQ(b1.FullName(), b2.FullName());
60     EXPECT_EQ(b1.HoleNumber(), b2.HoleNumber());
61     EXPECT_EQ(b1.NumPasses(), b2.NumPasses());
62     EXPECT_EQ(b1.Sequence(), b2.Sequence());
63     EXPECT_EQ(b1.Qualities(), b2.Qualities());
64     EXPECT_EQ(b1.DeletionQV(), b2.DeletionQV());
65     EXPECT_EQ(b1.DeletionTag(), b2.DeletionTag());
66     EXPECT_EQ(b1.InsertionQV(), b2.InsertionQV());
67     EXPECT_EQ(b1.MergeQV(), b2.MergeQV());
68     EXPECT_EQ(b1.SubstitutionQV(), b2.SubstitutionQV());
69     EXPECT_EQ(b1.SubstitutionTag(), b2.SubstitutionTag());
70     EXPECT_EQ(b1.LabelQV(), b2.LabelQV());
71     EXPECT_EQ(b1.AltLabelQV(), b2.AltLabelQV());
72     EXPECT_EQ(b1.AltLabelTag(), b2.AltLabelTag());
73     EXPECT_EQ(b1.Pkmean(), b2.Pkmean());
74     EXPECT_EQ(b1.Pkmid(), b2.Pkmid());
75     EXPECT_EQ(b1.PulseCall(), b2.PulseCall());
76     EXPECT_EQ(b1.IPD(), b2.IPD());
77     EXPECT_EQ(b1.PulseWidth(), b2.PulseWidth());
78     EXPECT_EQ(b1.PrePulseFrames(), b2.PrePulseFrames());
79     EXPECT_EQ(b1.PulseCallWidth(), b2.PulseCallWidth());
80     EXPECT_EQ(b1.ReadGroup(), b2.ReadGroup());
81     EXPECT_EQ(b1.PulseMergeQV(), b2.PulseMergeQV());
82 }
83 
NumVirtualRecords(const std::string & primaryBamFn,const std::string & scrapsBamFn)84 static size_t NumVirtualRecords(const std::string& primaryBamFn, const std::string& scrapsBamFn)
85 {
86     ZmwReadStitcher stitcher(primaryBamFn, scrapsBamFn);
87     size_t count = 0;
88     while (stitcher.HasNext()) {
89         const auto record = stitcher.Next();
90         //        ()record;
91         ++count;
92     }
93     return count;
94 }
95 
96 }  // namespace ZmwReadStitcherTests
97 
TEST(ZmwReadStitching,FromBams_NoFilter)98 TEST(ZmwReadStitching, FromBams_NoFilter)
99 {
100     ZmwReadStitcher stitcher(PbbamTestsConfig::Data_Dir + "/polymerase/internal.subreads.bam",
101                              PbbamTestsConfig::Data_Dir + "/polymerase/internal.scraps.bam");
102     size_t count = 0;
103     while (stitcher.HasNext()) {
104         const auto record = stitcher.Next();
105         //        ()record;
106         ++count;
107     }
108     EXPECT_EQ(3, count);
109 }
110 
TEST(ZmwReadStitching,FromBams_Filtered)111 TEST(ZmwReadStitching, FromBams_Filtered)
112 {
113     PbiFilter filter{PbiZmwFilter{100000}};  // setup to match DataSet w/ filter
114     ZmwReadStitcher stitcher(PbbamTestsConfig::Data_Dir + "/polymerase/internal.subreads.bam",
115                              PbbamTestsConfig::Data_Dir + "/polymerase/internal.scraps.bam",
116                              filter);
117     size_t count = 0;
118     while (stitcher.HasNext()) {
119         const auto record = stitcher.Next();
120         EXPECT_EQ(100000, record.HoleNumber());
121         ++count;
122     }
123     EXPECT_EQ(1, count);
124 }
125 
TEST(ZmwReadStitching,FromDataSet_NoFilter)126 TEST(ZmwReadStitching, FromDataSet_NoFilter)
127 {
128     // dataset contains these resources (subreads/scraps + hqregion/scraps BAMs)
129     const std::string primaryFn1 =
130         PbbamTestsConfig::Data_Dir + "/polymerase/production.subreads.bam";
131     const std::string scrapsFn1 = PbbamTestsConfig::Data_Dir + "/polymerase/production.scraps.bam";
132     const std::string primaryFn2 =
133         PbbamTestsConfig::Data_Dir + "/polymerase/production_hq.hqregion.bam";
134     const std::string scrapsFn2 =
135         PbbamTestsConfig::Data_Dir + "/polymerase/production_hq.scraps.bam";
136     const size_t numExpectedRecords =
137         ZmwReadStitcherTests::NumVirtualRecords(primaryFn1, scrapsFn1) +
138         ZmwReadStitcherTests::NumVirtualRecords(primaryFn2, scrapsFn2);
139 
140     const std::string datasetFn =
141         PbbamTestsConfig::Data_Dir + "/polymerase/multiple_resources.subread.dataset.xml";
142 
143     DataSet ds{datasetFn};
144     ZmwReadStitcher stitcher{ds};
145     size_t numObservedRecords = 0;
146     while (stitcher.HasNext()) {
147         const auto record = stitcher.Next();
148         //        ()record;
149         ++numObservedRecords;
150     }
151     EXPECT_EQ(numExpectedRecords, numObservedRecords);
152 }
153 
TEST(ZmwReadStitching,FromDataSet_Filtered)154 TEST(ZmwReadStitching, FromDataSet_Filtered)
155 {
156     // dataset contains these resources (subreads/scraps + hqregion/scraps BAMs)
157     const std::string primaryFn1 =
158         PbbamTestsConfig::Data_Dir + "/polymerase/production.subreads.bam";
159     const std::string scrapsFn1 = PbbamTestsConfig::Data_Dir + "/polymerase/production.scraps.bam";
160     const std::string primaryFn2 = PbbamTestsConfig::Data_Dir + "/polymerase/internal.subreads.bam";
161     const std::string scrapsFn2 = PbbamTestsConfig::Data_Dir + "/polymerase/internal.scraps.bam";
162     const std::string primaryFn3 =
163         PbbamTestsConfig::Data_Dir + "/polymerase/production_hq.hqregion.bam";
164     const std::string scrapsFn3 =
165         PbbamTestsConfig::Data_Dir + "/polymerase/production_hq.scraps.bam";
166     const size_t totalRecords = ZmwReadStitcherTests::NumVirtualRecords(primaryFn1, scrapsFn1) +
167                                 ZmwReadStitcherTests::NumVirtualRecords(primaryFn2, scrapsFn2) +
168                                 ZmwReadStitcherTests::NumVirtualRecords(primaryFn3, scrapsFn3);
169     EXPECT_EQ(5, totalRecords);
170 
171     // our filter will remove the 2 "production" BAM pairs
172     // using a ZMW filter that only the "internal" pair should pass
173     const std::string datasetFn =
174         PbbamTestsConfig::Data_Dir + "/polymerase/filtered_resources.subread.dataset.xml";
175 
176     DataSet ds{datasetFn};
177     ZmwReadStitcher stitcher{ds};
178     size_t numObservedRecords = 0;
179     while (stitcher.HasNext()) {
180         const auto record = stitcher.Next();
181         //        ()record;
182         ++numObservedRecords;
183     }
184     EXPECT_EQ(1, numObservedRecords);
185 }
186 
TEST(ZmwReadStitching,FromDataSet_EmptyDataSet)187 TEST(ZmwReadStitching, FromDataSet_EmptyDataSet)
188 {
189     ZmwReadStitcher stitcher{DataSet{}};
190     EXPECT_FALSE(stitcher.HasNext());
191 }
192 
TEST(ZmwReadStitching,EmptyScrapsFile)193 TEST(ZmwReadStitching, EmptyScrapsFile)
194 {
195     const std::string primaryBamFn =
196         PbbamTestsConfig::Data_Dir + "/polymerase/scrapless.subreads.bam";
197     const std::string scrapsBamFn = PbbamTestsConfig::Data_Dir + "/polymerase/scrapless.scraps.bam";
198 
199     const BamFile primaryBam(primaryBamFn);
200     const BamFile scrapsBam(scrapsBamFn);
201     const PbiRawData primaryIdx(primaryBam.PacBioIndexFilename());
202     const PbiRawData scrapsIdx(scrapsBam.PacBioIndexFilename());
203     EXPECT_EQ(3, primaryIdx.NumReads());
204     EXPECT_EQ(0, scrapsIdx.NumReads());
205 
206     int count = 0;
207     ZmwReadStitcher stitcher(primaryBamFn, scrapsBamFn);
208     while (stitcher.HasNext()) {
209         auto record = stitcher.Next();
210         //        ()record;
211         ++count;
212     }
213     EXPECT_EQ(3, count);
214 }
215 
TEST(ZmwReadStitching,VirtualRegions)216 TEST(ZmwReadStitching, VirtualRegions)
217 {
218     // Create virtual polymerase read
219     ZmwReadStitcher stitcher(PbbamTestsConfig::Data_Dir + "/polymerase/internal.subreads.bam",
220                              PbbamTestsConfig::Data_Dir + "/polymerase/internal.scraps.bam");
221     auto virtualRecord = stitcher.Next();
222 
223     auto regionMap = virtualRecord.VirtualRegionsMap();
224     auto adapter = virtualRecord.VirtualRegionsTable(VirtualRegionType::ADAPTER);
225 
226     // Compare different accessors to same source
227     EXPECT_EQ(regionMap[VirtualRegionType::ADAPTER], adapter);
228 
229     // Compare to truth
230     EXPECT_EQ(3047, adapter[0].beginPos);
231     EXPECT_EQ(3095, adapter[0].endPos);
232     EXPECT_EQ(3650, adapter[1].beginPos);
233     EXPECT_EQ(3700, adapter[1].endPos);
234     EXPECT_EQ(4289, adapter[2].beginPos);
235     EXPECT_EQ(4335, adapter[2].endPos);
236     EXPECT_EQ(4888, adapter[3].beginPos);
237     EXPECT_EQ(4939, adapter[3].endPos);
238     EXPECT_EQ(5498, adapter[4].beginPos);
239     EXPECT_EQ(5546, adapter[4].endPos);
240     EXPECT_EQ(6116, adapter[5].beginPos);
241     EXPECT_EQ(6173, adapter[5].endPos);
242     EXPECT_EQ(6740, adapter[6].beginPos);
243     EXPECT_EQ(6790, adapter[6].endPos);
244 
245     auto barcode = virtualRecord.VirtualRegionsTable(VirtualRegionType::BARCODE);
246     EXPECT_EQ(regionMap[VirtualRegionType::BARCODE], barcode);
247     EXPECT_EQ(3025, barcode[0].beginPos);
248     EXPECT_EQ(3047, barcode[0].endPos);
249     EXPECT_EQ(3095, barcode[1].beginPos);
250     EXPECT_EQ(3116, barcode[1].endPos);
251     EXPECT_EQ(3628, barcode[2].beginPos);
252     EXPECT_EQ(3650, barcode[2].endPos);
253     EXPECT_EQ(3700, barcode[3].beginPos);
254     EXPECT_EQ(3722, barcode[3].endPos);
255     EXPECT_EQ(4267, barcode[4].beginPos);
256     EXPECT_EQ(4289, barcode[4].endPos);
257     EXPECT_EQ(4335, barcode[5].beginPos);
258     EXPECT_EQ(4356, barcode[5].endPos);
259     EXPECT_EQ(4864, barcode[6].beginPos);
260     EXPECT_EQ(4888, barcode[6].endPos);
261     EXPECT_EQ(4939, barcode[7].beginPos);
262     EXPECT_EQ(4960, barcode[7].endPos);
263     EXPECT_EQ(5477, barcode[8].beginPos);
264     EXPECT_EQ(5498, barcode[8].endPos);
265     EXPECT_EQ(5546, barcode[9].beginPos);
266     EXPECT_EQ(5571, barcode[9].endPos);
267     EXPECT_EQ(6087, barcode[10].beginPos);
268     EXPECT_EQ(6116, barcode[10].endPos);
269     EXPECT_EQ(6173, barcode[11].beginPos);
270     EXPECT_EQ(6199, barcode[11].endPos);
271     EXPECT_EQ(6719, barcode[12].beginPos);
272     EXPECT_EQ(6740, barcode[12].endPos);
273     EXPECT_EQ(6790, barcode[13].beginPos);
274     EXPECT_EQ(6812, barcode[13].endPos);
275 
276     auto lqregion = virtualRecord.VirtualRegionsTable(VirtualRegionType::LQREGION);
277     EXPECT_EQ(regionMap[VirtualRegionType::LQREGION], lqregion);
278     EXPECT_EQ(0, lqregion[0].beginPos);
279     EXPECT_EQ(2659, lqregion[0].endPos);
280     EXPECT_EQ(7034, lqregion[1].beginPos);
281     EXPECT_EQ(7035, lqregion[1].endPos);
282 
283     auto hqregion = virtualRecord.VirtualRegionsTable(VirtualRegionType::HQREGION);
284     EXPECT_EQ(regionMap[VirtualRegionType::HQREGION], hqregion);
285     EXPECT_EQ(2659, hqregion[0].beginPos);
286     EXPECT_EQ(7034, hqregion[0].endPos);
287 }
288 
TEST(ZmwReadStitching,InternalSubreadsToOriginal)289 TEST(ZmwReadStitching, InternalSubreadsToOriginal)
290 {
291     // Create virtual polymerase read
292     ZmwReadStitcher stitcher(PbbamTestsConfig::Data_Dir + "/polymerase/internal.subreads.bam",
293                              PbbamTestsConfig::Data_Dir + "/polymerase/internal.scraps.bam");
294     EXPECT_TRUE(stitcher.HasNext());
295     auto virtualRecord = stitcher.Next();
296 
297     // Read original polymerase read
298     BamFile polyBam(PbbamTestsConfig::Data_Dir + "/polymerase/internal.polymerase.bam");
299     EntireFileQuery polyQuery(polyBam);
300     auto begin = polyQuery.begin();
301     auto end = polyQuery.end();
302     EXPECT_TRUE(begin != end);
303     auto polyRecord = *begin;
304 
305     // check
306     ZmwReadStitcherTests::Compare(polyRecord, virtualRecord);
307 }
308 
TEST(ZmwReadStitching,InternalHQToOriginal)309 TEST(ZmwReadStitching, InternalHQToOriginal)
310 {
311     // Create virtual polymerase read
312     ZmwReadStitcher stitcher(PbbamTestsConfig::Data_Dir + "/polymerase/internal.hqregions.bam",
313                              PbbamTestsConfig::Data_Dir + "/polymerase/internal.lqregions.bam");
314     EXPECT_TRUE(stitcher.HasNext());
315     auto virtualRecord = stitcher.Next();
316 
317     // Read original polymerase read
318     BamFile polyBam(PbbamTestsConfig::Data_Dir + "/polymerase/internal.polymerase.bam");
319     EntireFileQuery polyQuery(polyBam);
320     auto begin = polyQuery.begin();
321     auto end = polyQuery.end();
322     EXPECT_TRUE(begin != end);
323     auto polyRecord = *begin;
324 
325     // check
326     ZmwReadStitcherTests::Compare(polyRecord, virtualRecord);
327 }
328 
TEST(ZmwReadStitching,ProductionSubreadsToOriginal)329 TEST(ZmwReadStitching, ProductionSubreadsToOriginal)
330 {
331     // Create virtual polymerase read
332     ZmwReadStitcher stitcher(PbbamTestsConfig::Data_Dir + "/polymerase/production.subreads.bam",
333                              PbbamTestsConfig::Data_Dir + "/polymerase/production.scraps.bam");
334 
335     EXPECT_TRUE(stitcher.HasNext());
336     auto virtualRecord = stitcher.Next();
337     EXPECT_FALSE(stitcher.HasNext());
338 
339     // Read original polymerase read
340     BamFile polyBam(PbbamTestsConfig::Data_Dir + "/polymerase/production.polymerase.bam");
341     EntireFileQuery polyQuery(polyBam);
342 
343     auto begin = polyQuery.begin();
344     auto end = polyQuery.end();
345     EXPECT_TRUE(begin != end);
346     auto polyRecord = *begin;
347 
348     EXPECT_EQ(polyRecord.FullName(), virtualRecord.FullName());
349     EXPECT_EQ(polyRecord.HoleNumber(), virtualRecord.HoleNumber());
350     EXPECT_FLOAT_EQ(polyRecord.ReadAccuracy(), virtualRecord.ReadAccuracy());
351     EXPECT_EQ(polyRecord.NumPasses(), virtualRecord.NumPasses());
352     EXPECT_EQ(polyRecord.Sequence(), virtualRecord.Sequence());
353     EXPECT_EQ(polyRecord.Qualities(), virtualRecord.Qualities());
354     EXPECT_EQ(polyRecord.DeletionQV(), virtualRecord.DeletionQV());
355     EXPECT_EQ(polyRecord.DeletionTag(), virtualRecord.DeletionTag());
356     EXPECT_EQ(polyRecord.InsertionQV(), virtualRecord.InsertionQV());
357     EXPECT_EQ(polyRecord.MergeQV(), virtualRecord.MergeQV());
358     EXPECT_EQ(polyRecord.SubstitutionQV(), virtualRecord.SubstitutionQV());
359     EXPECT_EQ(polyRecord.SubstitutionTag(), virtualRecord.SubstitutionTag());
360     EXPECT_EQ(polyRecord.IPD(), virtualRecord.IPDV1Frames());
361     EXPECT_EQ(polyRecord.ReadGroup(), virtualRecord.ReadGroup());
362 }
363 
TEST(ZmwReadStitching,ProductionHQToOriginal)364 TEST(ZmwReadStitching, ProductionHQToOriginal)
365 {
366     // Create virtual polymerase read
367     ZmwReadStitcher stitcher(PbbamTestsConfig::Data_Dir + "/polymerase/production_hq.hqregion.bam",
368                              PbbamTestsConfig::Data_Dir + "/polymerase/production_hq.scraps.bam");
369     EXPECT_TRUE(stitcher.HasNext());
370     auto virtualRecord = stitcher.Next();
371     EXPECT_FALSE(stitcher.HasNext());
372 
373     // Read original polymerase read
374     BamFile polyBam(PbbamTestsConfig::Data_Dir + "/polymerase/production.polymerase.bam");
375     EntireFileQuery polyQuery(polyBam);
376 
377     auto begin = polyQuery.begin();
378     auto end = polyQuery.end();
379     EXPECT_TRUE(begin != end);
380     auto polyRecord = *begin;
381 
382     EXPECT_FALSE(polyRecord.HasPulseCall());
383     EXPECT_FALSE(virtualRecord.HasPulseCall());
384     EXPECT_EQ(polyRecord.FullName(), virtualRecord.FullName());
385     EXPECT_EQ(polyRecord.HoleNumber(), virtualRecord.HoleNumber());
386     EXPECT_EQ(polyRecord.ReadAccuracy(), virtualRecord.ReadAccuracy());
387     EXPECT_EQ(polyRecord.NumPasses(), virtualRecord.NumPasses());
388     EXPECT_EQ(polyRecord.Sequence(), virtualRecord.Sequence());
389     EXPECT_EQ(polyRecord.Qualities(), virtualRecord.Qualities());
390     EXPECT_EQ(polyRecord.DeletionQV(), virtualRecord.DeletionQV());
391     EXPECT_EQ(polyRecord.DeletionTag(), virtualRecord.DeletionTag());
392     EXPECT_EQ(polyRecord.InsertionQV(), virtualRecord.InsertionQV());
393     EXPECT_EQ(polyRecord.MergeQV(), virtualRecord.MergeQV());
394     EXPECT_EQ(polyRecord.SubstitutionQV(), virtualRecord.SubstitutionQV());
395     EXPECT_EQ(polyRecord.SubstitutionTag(), virtualRecord.SubstitutionTag());
396     EXPECT_EQ(polyRecord.IPD(), virtualRecord.IPDV1Frames());
397     EXPECT_EQ(polyRecord.ReadGroup(), virtualRecord.ReadGroup());
398 
399     EXPECT_TRUE(polyRecord.HasDeletionQV());
400     EXPECT_TRUE(polyRecord.HasDeletionTag());
401     EXPECT_TRUE(polyRecord.HasInsertionQV());
402     EXPECT_TRUE(polyRecord.HasMergeQV());
403     EXPECT_TRUE(polyRecord.HasSubstitutionQV());
404     EXPECT_TRUE(polyRecord.HasSubstitutionTag());
405     EXPECT_TRUE(polyRecord.HasIPD());
406     EXPECT_FALSE(polyRecord.HasLabelQV());
407     EXPECT_FALSE(polyRecord.HasAltLabelQV());
408     EXPECT_FALSE(polyRecord.HasAltLabelTag());
409     EXPECT_FALSE(polyRecord.HasPkmean());
410     EXPECT_FALSE(polyRecord.HasPkmid());
411     EXPECT_FALSE(polyRecord.HasPulseCall());
412     EXPECT_FALSE(polyRecord.HasPulseWidth());
413     EXPECT_FALSE(polyRecord.HasPrePulseFrames());
414     EXPECT_FALSE(polyRecord.HasPulseCallWidth());
415 
416     EXPECT_TRUE(virtualRecord.HasDeletionQV());
417     EXPECT_TRUE(virtualRecord.HasDeletionTag());
418     EXPECT_TRUE(virtualRecord.HasInsertionQV());
419     EXPECT_TRUE(virtualRecord.HasMergeQV());
420     EXPECT_TRUE(virtualRecord.HasSubstitutionQV());
421     EXPECT_TRUE(virtualRecord.HasSubstitutionTag());
422     EXPECT_TRUE(virtualRecord.HasIPD());
423     EXPECT_FALSE(virtualRecord.HasLabelQV());
424     EXPECT_FALSE(virtualRecord.HasAltLabelQV());
425     EXPECT_FALSE(virtualRecord.HasAltLabelTag());
426     EXPECT_FALSE(virtualRecord.HasPkmean());
427     EXPECT_FALSE(virtualRecord.HasPkmid());
428     EXPECT_FALSE(virtualRecord.HasPulseCall());
429     EXPECT_FALSE(virtualRecord.HasPulseWidth());
430     EXPECT_FALSE(virtualRecord.HasPrePulseFrames());
431     EXPECT_FALSE(virtualRecord.HasPulseCallWidth());
432 }
433 
TEST(ZmwReadStitching,VirtualRecord_VirtualRegionsTable)434 TEST(ZmwReadStitching, VirtualRecord_VirtualRegionsTable)
435 {
436     ZmwReadStitcher stitcher(PbbamTestsConfig::Data_Dir + "/polymerase/production.subreads.bam",
437                              PbbamTestsConfig::Data_Dir + "/polymerase/production.scraps.bam");
438     EXPECT_TRUE(stitcher.HasNext());
439     const auto virtualRecord = stitcher.Next();
440 
441     const auto subreads = virtualRecord.VirtualRegionsTable(VirtualRegionType::SUBREAD);
442     const auto adapters = virtualRecord.VirtualRegionsTable(VirtualRegionType::ADAPTER);
443     const auto hqRegions = virtualRecord.VirtualRegionsTable(VirtualRegionType::HQREGION);
444     const auto lqRegions = virtualRecord.VirtualRegionsTable(VirtualRegionType::LQREGION);
445     const auto barcodes = virtualRecord.VirtualRegionsTable(VirtualRegionType::BARCODE);
446     const auto filtered = virtualRecord.VirtualRegionsTable(VirtualRegionType::FILTERED);
447 
448     EXPECT_FALSE(subreads.empty());
449     EXPECT_FALSE(adapters.empty());
450     EXPECT_FALSE(hqRegions.empty());
451     EXPECT_FALSE(lqRegions.empty());
452     EXPECT_FALSE(barcodes.empty());
453     EXPECT_TRUE(filtered.empty());  // this type not present in this data
454 }
455 
TEST(ZmwReadStitching,LegacyTypedefsOk)456 TEST(ZmwReadStitching, LegacyTypedefsOk)
457 {
458     {
459         VirtualPolymeraseReader reader(
460             PbbamTestsConfig::Data_Dir + "/polymerase/internal.subreads.bam",
461             PbbamTestsConfig::Data_Dir + "/polymerase/internal.scraps.bam");
462         size_t count = 0;
463         while (reader.HasNext()) {
464             const auto record = reader.Next();
465             //            ()record;
466             ++count;
467         }
468         EXPECT_EQ(3, count);
469     }
470 
471     {
472         VirtualPolymeraseCompositeReader reader{DataSet{}};
473         EXPECT_FALSE(reader.HasNext());
474     }
475 }
476