htsworkflow/pipelines/test/simulate_runfolder.py

   1 """
   2 Create simulated solexa/illumina runfolders for testing
   3 """
   4 import gzip
   5 import os
   6 import shutil
   7
   8 TEST_CODE_DIR = os.path.split(__file__)[0]
   9 TESTDATA_DIR = os.path.join(TEST_CODE_DIR, 'testdata')
  10 LANE_LIST = range(1,9)
  11 TILE_LIST = range(1,101)
  12 HISEQ_TILE_LIST = [1101, 1102, 1103, 1104, 1105, 1106, 1107, 1108,
  13                    1201, 1202, 1203, 1204, 1205, 1206, 1207, 1208,
  14                    2101, 2102, 2103, 2104, 2105, 2106, 2107, 2108,
  15                    2201, 2202, 2203, 2204, 2205, 2206, 2207, 2208,]
  16
  17 def make_firecrest_dir(data_dir, version="1.9.2", start=1, stop=37):
  18     firecrest_dir = os.path.join(data_dir,
  19                                  'C%d-%d_Firecrest%s_12-04-2008_diane' % (start, stop, version)
  20                                  )
  21     os.mkdir(firecrest_dir)
  22     return firecrest_dir
  23
  24 def make_ipar_dir(data_dir, version='1.01'):
  25     """
  26     Construct an artificial ipar parameter file and directory
  27     """
  28     ipar1_01_file = os.path.join(TESTDATA_DIR, 'IPAR1.01.params')
  29     shutil.copy(ipar1_01_file, os.path.join(data_dir, '.params'))
  30
  31     ipar_dir = os.path.join(data_dir, 'IPAR_%s' % (version,))
  32     if not os.path.exists(ipar_dir):
  33       os.mkdir(ipar_dir)
  34     return ipar_dir
  35
  36 def make_flowcell_id(runfolder_dir, flowcell_id=None):
  37     if flowcell_id is None:
  38         flowcell_id = '207BTAAXY'
  39
  40     config = """<?xml version="1.0"?>
  41 <FlowcellId>
  42   <Text>%s</Text>
  43 </FlowcellId>""" % (flowcell_id,)
  44     config_dir = os.path.join(runfolder_dir, 'Config')
  45
  46     if not os.path.exists(config_dir):
  47         os.mkdir(config_dir)
  48     pathname = os.path.join(config_dir, 'FlowcellId.xml')
  49     f = open(pathname,'w')
  50     f.write(config)
  51     f.close()
  52
  53 def make_bustard_config132(image_dir):
  54     source = os.path.join(TESTDATA_DIR, 'bustard-config132.xml')
  55     destination = os.path.join(image_dir, 'config.xml')
  56     shutil.copy(source, destination)
  57
  58 def make_aligned_config_1_12(aligned_dir):
  59     """This is rouglhly equivalent to the old gerald file"""
  60     source = os.path.join(TESTDATA_DIR, '1_12', 'aligned_config_1_12.xml')
  61     destination = os.path.join(aligned_dir, 'config.xml')
  62     shutil.copy(source, destination)
  63
  64 def make_unaligned_config_1_12(unaligned_dir):
  65     demultiplex_pairs = [ # (src,
  66       # dest),
  67         (os.path.join(TESTDATA_DIR, '1_12', 'demultiplex_1.12.4.2.xml'),
  68          os.path.join(unaligned_dir, 'DemultiplexConfig.xml')),
  69         (os.path.join(TESTDATA_DIR, '1_12',
  70                       'demultiplexed_bustard_1.12.4.2.xml'),
  71          os.path.join(unaligned_dir, 'DemultiplexedBustardConfig.xml')),
  72         (os.path.join(TESTDATA_DIR, '1_12',
  73                       'demultiplexed_summary_1.12.4.2.xml'),
  74          os.path.join(unaligned_dir, 'DemultiplexedBustardSummary.xml')),
  75     ]
  76     for src, dest in demultiplex_pairs:
  77         shutil.copy(src, dest)
  78
  79 def make_rta_intensities_1460(data_dir, version='1.4.6.0'):
  80     """
  81     Construct an artificial RTA Intensities parameter file and directory
  82     """
  83     intensities_dir = os.path.join(data_dir, 'Intensities')
  84     if not os.path.exists(intensities_dir):
  85       os.mkdir(intensities_dir)
  86
  87     param_file = os.path.join(TESTDATA_DIR, 'rta_intensities_config.xml')
  88     shutil.copy(param_file, os.path.join(intensities_dir, 'config.xml'))
  89
  90     return intensities_dir
  91
  92 def make_rta_basecalls_1460(intensities_dir):
  93     """
  94     Construct an artificial RTA Intensities parameter file and directory
  95     """
  96     basecalls_dir = os.path.join(intensities_dir, 'BaseCalls')
  97     if not os.path.exists(basecalls_dir):
  98       os.mkdir(basecalls_dir)
  99
 100     param_file = os.path.join(TESTDATA_DIR, 'rta_basecalls_config.xml')
 101     shutil.copy(param_file, os.path.join(basecalls_dir, 'config.xml'))
 102
 103     return basecalls_dir
 104
 105 def make_rta_intensities_1870(data_dir, version='1.8.70.0'):
 106     """
 107     Construct an artificial RTA Intensities parameter file and directory
 108     """
 109     intensities_dir = os.path.join(data_dir, 'Intensities')
 110     if not os.path.exists(intensities_dir):
 111       os.mkdir(intensities_dir)
 112
 113     param_file = os.path.join(TESTDATA_DIR, 'rta_intensities_config_1870.xml')
 114     shutil.copy(param_file, os.path.join(intensities_dir, 'config.xml'))
 115
 116     return intensities_dir
 117
 118 def make_rta_intensities_1_10(data_dir, version='1.10.36.0'):
 119     """
 120     Construct an artificial RTA Intensities parameter file and directory
 121     """
 122     intensities_dir = os.path.join(data_dir, 'Intensities')
 123     if not os.path.exists(intensities_dir):
 124       os.mkdir(intensities_dir)
 125
 126     param_file = os.path.join(TESTDATA_DIR, 'rta_intensities_config_1.10.xml')
 127     shutil.copy(param_file, os.path.join(intensities_dir, 'config.xml'))
 128
 129     return intensities_dir
 130
 131 def make_rta_intensities_1_12(data_dir, version='1.12.4.2'):
 132     """
 133     Construct an artificial RTA Intensities parameter file and directory
 134     """
 135     intensities_dir = os.path.join(data_dir, 'Intensities')
 136     if not os.path.exists(intensities_dir):
 137       os.mkdir(intensities_dir)
 138
 139     param_file = os.path.join(TESTDATA_DIR, '1_12',
 140                               'rta_intensities_config_1.12.4.2.xml')
 141     shutil.copy(param_file, os.path.join(intensities_dir, 'RTAConfig.xml'))
 142
 143     return intensities_dir
 144
 145 def make_rta_basecalls_1870(intensities_dir):
 146     """
 147     Construct an artificial RTA Intensities parameter file and directory
 148     """
 149     basecalls_dir = os.path.join(intensities_dir, 'BaseCalls')
 150     if not os.path.exists(basecalls_dir):
 151       os.mkdir(basecalls_dir)
 152
 153     param_file = os.path.join(TESTDATA_DIR, 'rta_basecalls_config_1870.xml')
 154     shutil.copy(param_file, os.path.join(basecalls_dir, 'config.xml'))
 155
 156     return basecalls_dir
 157
 158 def make_rta_basecalls_1_10(intensities_dir):
 159     """
 160     Construct an artificial RTA Intensities parameter file and directory
 161     """
 162     basecalls_dir = os.path.join(intensities_dir, 'BaseCalls')
 163     if not os.path.exists(basecalls_dir):
 164         os.mkdir(basecalls_dir)
 165
 166     make_qseqs(basecalls_dir, basecall_info=ABXX_BASE_CALL_INFO)
 167     param_file = os.path.join(TESTDATA_DIR, 'rta_basecalls_config_1.10.xml')
 168     shutil.copy(param_file, os.path.join(basecalls_dir, 'config.xml'))
 169
 170     return basecalls_dir
 171
 172 def make_rta_basecalls_1_12(intensities_dir):
 173     """
 174     Construct an artificial RTA Intensities parameter file and directory
 175     """
 176     basecalls_dir = os.path.join(intensities_dir, 'BaseCalls')
 177     if not os.path.exists(basecalls_dir):
 178         os.mkdir(basecalls_dir)
 179
 180     make_qseqs(basecalls_dir, basecall_info=ABXX_BASE_CALL_INFO)
 181     param_file = os.path.join(TESTDATA_DIR, '1_12',
 182                               'rta_basecalls_config_1.12.4.2.xml')
 183     shutil.copy(param_file, os.path.join(basecalls_dir, 'config.xml'))
 184
 185     return basecalls_dir
 186
 187
 188 def make_qseqs(bustard_dir, basecall_info=None):
 189     """
 190     Fill gerald directory with qseq files
 191     """
 192     if basecall_info is None:
 193         qseq_file = '42BRJAAXX_8_1_0039_qseq.txt'
 194         tile_list = TILE_LIST
 195         summary_file = '42BRJAAXX_BustardSummary.xml'
 196     else:
 197         qseq_file = basecall_info.qseq_file
 198         tile_list = basecall_info.tile_list
 199         summary_file = basecall_info.basecall_summary
 200
 201     # 42BRJ 8 1 0039 happened to be a better than usual tile, in that there
 202     # was actually sequence at the start
 203     source = os.path.join(TESTDATA_DIR, qseq_file)
 204     destdir = bustard_dir
 205     if not os.path.isdir(destdir):
 206         os.mkdir(destdir)
 207
 208     for lane in LANE_LIST:
 209         for tile in tile_list:
 210             destination = os.path.join(bustard_dir, 's_%d_1_%04d_qseq.txt' % (lane, tile))
 211             shutil.copy(source, destination)
 212
 213     make_matrix_dir(bustard_dir)
 214     make_phasing_dir(bustard_dir)
 215
 216     summary_source = os.path.join(TESTDATA_DIR, summary_file)
 217     summary_dest = os.path.join(bustard_dir, 'BustardSummary.xml')
 218     shutil.copy(summary_source, summary_dest)
 219
 220     return destdir
 221
 222 def make_scores(gerald_dir, in_temp=True):
 223     """
 224     Fill gerald directory with score temp files
 225     will create the directory if it doesn't exist.
 226     """
 227     source = os.path.join(TESTDATA_DIR, 's_1_0001_score.txt')
 228     destdir = gerald_dir
 229     if in_temp:
 230         destdir = os.path.join(destdir, 'Temp')
 231     if not os.path.isdir(destdir):
 232         os.mkdir(destdir)
 233
 234     for lane in LANE_LIST:
 235         for tile in TILE_LIST:
 236             destination = os.path.join(destdir, 's_%d_%04d_score.txt' % (lane, tile))
 237             shutil.copy(source, destination)
 238
 239     return destdir
 240
 241 def make_matrix_dir(bustard_dir):
 242     """
 243     Create several matrix files in <bustard_dir>/Matrix/
 244
 245     from pipeline 1.4
 246     """
 247     destdir = os.path.join(bustard_dir, 'Matrix')
 248     if not os.path.isdir(destdir):
 249         os.mkdir(destdir)
 250
 251     source = os.path.join(TESTDATA_DIR, '42BRJAAXX_8_02_matrix.txt')
 252     for lane in LANE_LIST:
 253         destination = os.path.join(destdir, 's_%d_02_matrix.txt' % ( lane, ))
 254         shutil.copy(source, destination)
 255
 256 def make_matrix(matrix_filename):
 257     contents = """# Auto-generated frequency response matrix
 258 > A
 259 > C
 260 > G
 261 > T
 262 0.77 0.15 -0.04 -0.04
 263 0.76 1.02 -0.05 -0.06
 264 -0.10 -0.10 1.17 -0.03
 265 -0.13 -0.12 0.80 1.27
 266 """
 267     f = open(matrix_filename, 'w')
 268     f.write(contents)
 269     f.close()
 270
 271 def make_matrix_dir_rta160(bustard_dir):
 272     """
 273     Create several matrix files in <bustard_dir>/Matrix/
 274     """
 275     destdir = os.path.join(bustard_dir, 'Matrix')
 276     if not os.path.isdir(destdir):
 277         os.mkdir(destdir)
 278
 279     source = os.path.join(TESTDATA_DIR, '61MMFAAXX_4_1_matrix.txt')
 280     lane_fragments = [ "_%d" % (l,) for l in LANE_LIST]
 281     for fragment in lane_fragments:
 282         destination = os.path.join(destdir, 's%s_1_matrix.txt' % ( fragment, ))
 283         shutil.copy(source, destination)
 284
 285 def make_matrix_dir_rta_1_10(bustard_dir):
 286     make_matrix_dir_rta160(bustard_dir)
 287
 288 def make_matrix_dir_rta_1_12(bustard_dir):
 289     make_matrix_dir_rta160(bustard_dir)
 290
 291 def make_phasing_dir(bustard_dir):
 292     """
 293     Create several phasing files in <bustard_dir>/Phasing/
 294
 295     from pipeline 1.4
 296     """
 297     destdir = os.path.join(bustard_dir, 'Phasing')
 298     if not os.path.isdir(destdir):
 299         os.mkdir(destdir)
 300
 301     source = os.path.join(TESTDATA_DIR, '42BRJAAXX_8_01_phasing.xml')
 302     for lane in LANE_LIST:
 303         destination = os.path.join(destdir, 's_%d_01_phasing.xml' % ( lane, ))
 304         shutil.copy(source, destination)
 305
 306 def make_phasing_params(bustard_dir):
 307     for lane in LANE_LIST:
 308         pathname = os.path.join(bustard_dir, 'params%d.xml' % (lane))
 309         f = open(pathname, 'w')
 310         f.write("""<Parameters>
 311   <Phasing>0.009900</Phasing>
 312   <Prephasing>0.003500</Prephasing>
 313 </Parameters>
 314 """)
 315         f.close()
 316
 317 def make_gerald_config_026(gerald_dir):
 318     source = os.path.join(TESTDATA_DIR, 'gerald_config_0.2.6.xml')
 319     destination = os.path.join(gerald_dir, 'config.xml')
 320     shutil.copy(source, destination)
 321
 322 def make_gerald_config_100(gerald_dir):
 323     source = os.path.join(TESTDATA_DIR, 'gerald_config_1.0.xml')
 324     destination = os.path.join(gerald_dir, 'config.xml')
 325     shutil.copy(source, destination)
 326
 327 def make_gerald_config_1_7(gerald_dir):
 328     """CASAVA 1.7 gerald config"""
 329     source = os.path.join(TESTDATA_DIR, 'gerald_config_1.7.xml')
 330     destination = os.path.join(gerald_dir, 'config.xml')
 331     shutil.copy(source, destination)
 332
 333 def make_summary_htm_100(gerald_dir):
 334     source = os.path.join(TESTDATA_DIR, 'Summary-pipeline100.htm')
 335     destination = os.path.join(gerald_dir, 'Summary.htm')
 336     shutil.copy(source, destination)
 337
 338 def make_summary_htm_110(gerald_dir):
 339     source = os.path.join(TESTDATA_DIR, 'Summary-pipeline110.htm')
 340     destination = os.path.join(gerald_dir, 'Summary.htm')
 341     shutil.copy(source, destination)
 342
 343 def make_summary_paired_htm(gerald_dir):
 344     source = os.path.join(TESTDATA_DIR, 'Summary-paired-pipeline110.htm')
 345     destination = os.path.join(gerald_dir, 'Summary.htm')
 346     shutil.copy(source, destination)
 347
 348 def make_summary_ipar130_htm(gerald_dir):
 349     source = os.path.join(TESTDATA_DIR, 'Summary-ipar130.htm')
 350     destination = os.path.join(gerald_dir, 'Summary.htm')
 351     shutil.copy(source, destination)
 352
 353 def make_summary_rta160_xml(gerald_dir):
 354     source = os.path.join(TESTDATA_DIR, 'Summary-rta160.xml')
 355     destination = os.path.join(gerald_dir, 'Summary.xml')
 356     shutil.copy(source, destination)
 357
 358
 359 def make_summary_casava1_7_xml(gerald_dir):
 360     source = os.path.join(TESTDATA_DIR, 'Summary-casava1.7.xml')
 361     destination = os.path.join(gerald_dir, 'Summary.xml')
 362     shutil.copy(source, destination)
 363
 364 def make_status_rta1_12(datadir):
 365     sourcedir = os.path.join(TESTDATA_DIR, '1_12')
 366     status_htm = os.path.join(sourcedir, 'Status.htm')
 367     destination = os.path.join(datadir, 'Status.htm')
 368     shutil.copy(status_htm, destination)
 369
 370     status_dir = os.path.join(datadir, 'Status_Files')
 371     status_source_dir = os.path.join(sourcedir, 'Status_Files')
 372     shutil.copytree(status_source_dir, status_dir)
 373
 374     report_source_dir = os.path.join(sourcedir, 'reports')
 375     report_dir = os.path.join(datadir, 'reports')
 376     shutil.copytree(report_source_dir, report_dir)
 377
 378 def make_eland_results(gerald_dir):
 379     eland_result = """>HWI-EAS229_24_207BTAAXX:1:7:599:759    ACATAGNCACAGACATAAACATAGACATAGAC U0      1       1       3       chrUextra.fa    28189829        R       D.
 380 >HWI-EAS229_24_207BTAAXX:1:7:205:842    AAACAANNCTCCCAAACACGTAAACTGGAAAA  U1      0       1       0       chr2L.fa        8796855 R       DD      24T
 381 >HWI-EAS229_24_207BTAAXX:1:7:776:582    AGCTCANCCGATCGAAAACCTCNCCAAGCAAT        NM      0       0       0
 382 >HWI-EAS229_24_207BTAAXX:1:7:205:842    AAACAANNCTCCCAAACACGTAAACTGGAAAA        U1      0       1       0       Lambda.fa        8796855 R       DD      24T
 383 """
 384     for i in LANE_LIST:
 385         pathname = os.path.join(gerald_dir,
 386                                 's_%d_eland_result.txt' % (i,))
 387         f = open(pathname, 'w')
 388         f.write(eland_result)
 389         f.close()
 390
 391 def make_eland_multi(gerald_dir, paired=False, lane_list=LANE_LIST):
 392     eland_multi = [""">HWI-EAS229_60_30DP9AAXX:1:1:1221:788   AAGATATCTACGACGTGGTATGGCGGTGTCTGGTCGT      NM
 393 >HWI-EAS229_60_30DP9AAXX:1:1:931:747    AAAAAAGCAAATTTCATTCACATGTTCTGTGTTCATA   1:0:2   chr5.fa:55269838R0
 394 >HWI-EAS229_60_30DP9AAXX:1:1:1121:379   AGAAGAGACATTAAGAGTTCCTGAAATTTATATCTGG   2:1:0   chr16.fa:46189180R1,chr7.fa:122968519R0,chr8.fa:48197174F0
 395 >HWI-EAS229_60_30DP9AAXX:1:1:892:1155   ACATTCTCCTTTCCTTCTGAAGTTTTTACGATTCTTT   0:9:10  chr10.fa:114298201F1,chr12.fa:8125072F1,19500297F2,42341293R2,chr13.fa:27688155R2,95069772R1,chr15.fa:51016475F2,chr16.fa:27052155F2,chr1.fa:192426217R2,chr21.fa:23685310R2,chr2.fa:106680068F1,chr3.fa:185226695F2,chr4.fa:106626808R2,chr5.fa:14704894F1,43530779F1,126543189F2,chr6.fa:74284101F1,chr7.fa:22516603F1,chr9.fa:134886204R
 396 >HWI-EAS229_60_30DP9AAXX:1:1:931:747    AAAAAAGCAAATTTCATTCACATGTTCTGTGTTCATA   1:0:0   spike.fa/sample1:55269838R0
 397 >HWI-EAS229_60_30DP9AAXX:1:1:931:747    AAAAAAGCAAATTTCATTCACATGTTCTGTGTTCATA   1:0:0   spike.fa/sample2:55269838R0
 398 """, """>HWI-EAS229_60_30DP9AAXX:1:1:1221:788   AAGATATCTACGACGTGGTATGGCGGTGTCTGGTCGT      NM
 399 >HWI-EAS229_60_30DP9AAXX:1:1:1221:788   NNNNNNNNNNNNNNGTGGTATGGCGGTGTCTGGTCGT     QC
 400 >HWI-EAS229_60_30DP9AAXX:1:1:931:747    AAAAAAGCAAATTTCATTCACATGTTCTGTGTTCATA   1:0:2   chr5.fa:55269838R0
 401 >HWI-EAS229_60_30DP9AAXX:1:1:1121:379   AGAAGAGACATTAAGAGTTCCTGAAATTTATATCTGG   2:1:0   chr16.fa:46189180R1,chr7.fa:122968519R0,chr8.fa:48197174F0,chr7.fa:22516603F1,chr9.fa:134886204R
 402 >HWI-EAS229_60_30DP9AAXX:1:1:892:1155   ACATTCTCCTTTCCTTCTGAAGTTTTTACGATTCTTT   0:9:10  chr10.fa:114298201F1,chr12.fa:8125072F1,19500297F2,42341293R2,chr13.fa:27688155R2,95069772R1,chr15.fa:51016475F2,chr16.fa:27052155F2,chr1.fa:192426217R2,chr21.fa:23685310R2,chr2.fa:106680068F1,chr3.fa:185226695F2,chr4.fa:106626808R2,chr5.fa:14704894F1,43530779F1,126543189F2,chr6.fa:74284101F1
 403 >HWI-EAS229_60_30DP9AAXX:1:1:931:747    AAAAAAGCAAATTTCATTCACATGTTCTGTGTTCATA   1:0:0   spike.fa/sample1:55269838R0
 404 >HWI-EAS229_60_30DP9AAXX:1:1:931:747    AAAAAAGCAAATTTCATTCACATGTTCTGTGTTCATA   1:0:0   spike.fa/sample2:55269838R0
 405 """]
 406     if paired:
 407         for e in [1,2]:
 408             for i in lane_list:
 409                 pathname = os.path.join(gerald_dir,
 410                                         's_%d_%d_eland_multi.txt' % (i,e))
 411                 f = open(pathname, 'w')
 412                 f.write(eland_multi[e-1])
 413                 f.close()
 414     else:
 415         for i in lane_list:
 416             pathname = os.path.join(gerald_dir,
 417                                     's_%d_eland_multi.txt' % (i,))
 418             f = open(pathname, 'w')
 419             f.write(eland_multi[0])
 420             f.close()
 421
 422 def make_eland_export(gerald_dir, paired=False, lane_list=LANE_LIST):
 423     source = os.path.join(TESTDATA_DIR, 'casava_1.7_export.txt')
 424
 425     for i in lane_list:
 426         destination = os.path.join(gerald_dir,
 427                                    's_%d_export.txt' % (i,))
 428         shutil.copy(source, destination)
 429
 430
 431 def make_scarf(gerald_dir, lane_list=LANE_LIST):
 432     seq = """HWI-EAS229_92_30VNBAAXX:1:1:0:161:NCAATTACACGACGCTAGCCCTAAAGCTATTTCGAGG:E[aaaabb^a\a_^^a[S`ba_WZUXaaaaaaUKPER
 433 HWI-EAS229_92_30VNBAAXX:1:1:0:447:NAGATGCGCATTTGAAGTAGGAGCAAAAGATCAAGGT:EUabaab^baabaaaaaaaa^^Uaaaaa\aaaa__`a
 434 HWI-EAS229_92_30VNBAAXX:1:1:0:1210:NATAGCCTCTATAGAAGCCACTATTATTTTTTTCTTA:EUa`]`baaaaa^XQU^a`S``S_`J_aaaaaabb^V
 435 HWI-EAS229_92_30VNBAAXX:1:1:0:1867:NTGGAGCAGATATAAAAACAGATGGTGACGTTGAAGT:E[^UaaaUaba^aaa^aa^XV\baaLaLaaaaQVXV^
 436 HWI-EAS229_92_30VNBAAXX:1:1:0:1898:NAGCTCGTGTCGTGAGATGTTAGGTTAAGTCCTGCAA:EK_aaaaaaaaaaaUZaaZaXM[aaaXSM\aaZ]URE
 437 """
 438     for l in lane_list:
 439         pathname = os.path.join(gerald_dir, 's_%d_sequence.txt' %(l,))
 440         f = open(pathname,'w')
 441         f.write(seq)
 442         f.close()
 443
 444 def make_fastq(gerald_dir, lane_list=LANE_LIST):
 445     seq = """@HWI-EAS229:1:2:182:712#0/1
 446 AAAAAAAAAAAAAAAAAAAAANAAAAAAAAAAAAAAA
 447 +HWI-EAS229:1:2:182:712#0/1
 448 \\bab_bbaabbababbaaa]]D]bb_baabbab\baa
 449 @HWI-EAS229:1:2:198:621#0/1
 450 CCCCCCCCCCCCCCCCCCCCCNCCCCCCCCCCCCCCC
 451 +HWI-EAS229:1:2:198:621#0/1
 452 [aaaaaaa`_`aaaaaaa[`ZDZaaaaaaaaaaaaaa
 453 @HWI-EAS229:1:2:209:1321#0/1
 454 AAAAAAAAAAAAAAAAAAAAANAAAAAAAAAAAAAAA
 455 +HWI-EAS229:1:2:209:1321#0/1
 456 _bbbbbaaababaabbbbab]D]aaaaaaaaaaaaaa
 457 """
 458     for l in lane_list:
 459         pathname = os.path.join(gerald_dir, 's_%d_sequence.txt' %(l,))
 460         f = open(pathname,'w')
 461         f.write(seq)
 462         f.close()
 463
 464 UNALIGNED_READS = [1,2]
 465 UNALIGNED_SAMPLES = [ (1, UNALIGNED_READS, '11111', None, None),
 466                       (2, UNALIGNED_READS, '11112', None, None),
 467                       (3, UNALIGNED_READS, '11113', 1, 'ATCACG'),
 468                       (3, UNALIGNED_READS, '11113', 2, 'CGATGT'),
 469                       (3, UNALIGNED_READS, '11113', 3, 'TTAGGC'),
 470                       (4, UNALIGNED_READS, '11114', 6, 'GCCAAT'),
 471                       (5, UNALIGNED_READS, '11115', 1, 'ATCACG'),
 472                       (5, UNALIGNED_READS, '11116', 7, 'ACTTGA'),
 473                       (5, UNALIGNED_READS, '11117', 9, 'GATCAG'),
 474                       (6, UNALIGNED_READS, '11118', 1, 'ATCACG'),
 475                       (7, UNALIGNED_READS, '11119', 2, 'CGATGT'),
 476                       (8, UNALIGNED_READS, '11120', 3, 'TTAGGC'),
 477                       (1, UNALIGNED_READS, None, None, None),
 478                       (2, UNALIGNED_READS, None, None, None),
 479                       (3, UNALIGNED_READS, None, None, None),
 480                       (4, UNALIGNED_READS, None, None, None),
 481                       (5, UNALIGNED_READS, None, None, None)]
 482
 483
 484 def make_aligned_eland_export(aligned_dir, flowcell_id):
 485     summary_source = os.path.join(TESTDATA_DIR, 'sample_summary_1_12.htm')
 486     for lane, read, project_id, index_id, index_seq in UNALIGNED_SAMPLES:
 487         paths = DemultiplexedPaths(aligned_dir,
 488                                    flowcell_id,
 489                                    lane,
 490                                    project_id,
 491                                    index_id,
 492                                    index_seq)
 493         paths.make_sample_dirs()
 494         paths.make_summary_dirs()
 495         summary_dest = os.path.join(paths.summary_dir, 'Sample_Summary.htm')
 496         shutil.copy(summary_source, summary_dest)
 497
 498         body = get_unaligned_sample_export(lane, index_seq)
 499         for split in ['001','002']:
 500             for read in UNALIGNED_READS:
 501                 suffix = 'R{0}_{1}_export.txt.gz'.format(read, split)
 502                 pathname = paths.make_test_filename(suffix)
 503                 stream = gzip.open(pathname, 'w')
 504                 stream.write(body)
 505                 stream.close()
 506
 507
 508 def make_unaligned_fastqs_1_12(unaligned_dir, flowcell_id):
 509     """Create a default mix of unaligned sample files
 510     """
 511     for lane, read, name, index_id, index in UNALIGNED_SAMPLES:
 512         make_unaligned_fastq_sample_1_12(unaligned_dir,
 513                                          flowcell_id,
 514                                          lane,
 515                                          read,
 516                                          name,
 517                                          index_id,
 518                                          index)
 519
 520 def make_unaligned_fastq_sample_1_12(unaligned_dir,
 521                                      flowcell_id,
 522                                      lane,
 523                                      reads,
 524                                      project_id,
 525                                      index_id=None,
 526                                      index_seq=None):
 527
 528     paths = DemultiplexedPaths(unaligned_dir,
 529                                flowcell_id,
 530                                lane,
 531                                project_id,
 532                                index_id,
 533                                index_seq)
 534     paths.make_sample_dirs()
 535
 536     sample_seq = get_unaligned_sample_fastq_data(flowcell_id, lane, index_seq)
 537     for split in ['001','002']:
 538         for read in reads:
 539             suffix = 'R{0}_{1}.fastq.gz'.format(read, split)
 540             pathname = paths.make_test_filename(suffix)
 541             stream = gzip.open(pathname, 'w')
 542             stream.write(sample_seq)
 543             stream.close()
 544
 545     sheetname = os.path.join(paths.sample_dir, 'SampleSheet.csv')
 546     stream = open(sheetname, 'w')
 547     stream.write('FCID,Lane,SampleID,SampleRef,Index,Description,Control,Recipe,Operator,SampleProject'+os.linesep)
 548     template = '{flowcell},{lane},{id},mm9,{index},Sample #{id},N,PR_indexing,Operator,{sample_project}'+os.linesep
 549     stream.write(template.format(flowcell=flowcell_id,
 550                                  lane=lane,
 551                                  id=paths.sample_id,
 552                                  index=paths.index_seq,
 553                                  sample_project=paths.sample_project))
 554     stream.close()
 555
 556
 557 class DemultiplexedPaths(object):
 558     def __init__(self, basedir, flowcell_id, lane, project_id, index_id, index_seq):
 559         if lane not in LANE_LIST:
 560             raise ValueError("Invalid lane ID: {0}".format(lane))
 561         self.basedir = basedir
 562         self.flowcell_id = flowcell_id
 563         self.lane = lane
 564
 565         if project_id is None:
 566             # undetermined
 567             self.index_seq = ''
 568             self.sample_id = 'lane{0}'.format(lane)
 569             self.sample_project = 'Undetermined_indices'
 570             self.rootname = 'lane{lane}_Undetermined_L00{lane}_'.format(
 571                 lane=lane)
 572             self.project_dir = 'Undetermined_indices'
 573             self.sample_dir = 'Sample_lane{lane}'.format(lane=lane)
 574         elif index_seq is None:
 575             self.index_seq = ''
 576             self.sample_id = project_id
 577             self.sample_project = '{project_id}'.format(project_id=project_id)
 578             self.rootname = '{project_id}_NoIndex_L00{lane}_'.format(
 579                 project_id=project_id,
 580                 lane=lane)
 581             self.project_dir = 'Project_' + self.sample_project
 582             self.sample_dir = 'Sample_{project_id}'.format(
 583                 project_id=project_id)
 584         else:
 585             self.index_seq = index_seq
 586             self.sample_id = project_id
 587             self.sample_project = '{project_id}_Index{index_id}'.format(
 588                 project_id=project_id,
 589                 index_id=index_id)
 590             self.rootname = '{project_id}_{index}_L00{lane}_'.format(
 591                 project_id=project_id,
 592                 index=index_seq,
 593                 lane=lane)
 594             self.project_dir = 'Project_' + self.sample_project
 595             self.sample_dir = 'Sample_{project_id}'.format(
 596                 project_id=project_id)
 597
 598         self.project_dir = os.path.join(self.basedir, self.project_dir)
 599         self.sample_dir = os.path.join(self.project_dir, self.sample_dir)
 600         self.summary_dir = 'Summary_Stats_{0}'.format(self.flowcell_id)
 601         self.summary_dir = os.path.join(self.project_dir, self.summary_dir)
 602
 603
 604     def make_sample_dirs(self):
 605         if not os.path.isdir(self.project_dir):
 606             os.mkdir(self.project_dir)
 607         if not os.path.isdir(self.sample_dir):
 608             os.mkdir(self.sample_dir)
 609
 610     def make_summary_dirs(self):
 611         if not os.path.isdir(self.summary_dir):
 612             os.mkdir(self.summary_dir)
 613
 614     def make_test_filename(self, suffix):
 615         filename = self.rootname + suffix
 616         pathname = os.path.join(self.sample_dir, filename)
 617         return pathname
 618
 619     def dump(self):
 620         print ('index seq: {0}'.format(self.index_seq))
 621
 622         print ('project dir: {0}'.format(self.project_dir))
 623         print ('sample dir: {0}'.format(self.sample_dir))
 624         print ('rootname: {0}'.format(self.rootname))
 625         print ('path: {0}'.format(
 626             os.path.join(self.project_dir,
 627                          self.sample_dir,
 628                          self.rootname+'R1_001.fastq.gz')))
 629
 630
 631 def get_unaligned_sample_fastq_data(flowcell_id, lane, index_seq):
 632     seq = """@HWI-ST0787:101:{flowcell}:{lane}:1101:2416:3469 1:Y:0:{index}
 633 TCCTTCATTCCACCGGAGTCTGTGGAATTCTCGGGTGCCAAGGAACTCCA
 634 +
 635 CCCFFFFFHHHHHJJJJJJJJJIJJJJJJJJJJJJJJJJIIJJIIJJJJJ
 636 @HWI-ST0787:101:{flowcell}:{lane}:1101:2677:3293 1:Y:0:{index}
 637 TGGAAATCCATTGGGGTTTCCCCTGGAATTCTCGGGTGCCAAGGAACTCC
 638 +
 639 @CCFF3BDHHHHHIIIIIHHIIIDIIIGIIIEGIIIIIIIIIIIIIIIHH
 640 @HWI-ST0787:101:{flowcell}:{lane}:1101:2616:3297 1:Y:0:{index}
 641 TAATACTGCCGGGTAATGATGGCTGGAATTCTCGGGTGCCAAGGAACTCC
 642 +
 643 CCCFFFFFHHHHHCGHJJJJJJJJJJJJJJJJJIIJJJJJJJJJIHJJJI
 644 @HWI-ST0787:101:{flowcell}:{lane}:1101:2545:3319 1:N:0:{index}
 645 TCCTTCATTCCACCGGAGTCTGCTGGAATTCTCGGGTGCCAAGGAACTCC
 646 +
 647 CCCFFFFFHHHFHJGIGHIJHIIGHIGIGIGEHFIJJJIHIJHJIIJJIH
 648 """.format(flowcell=flowcell_id, lane=lane, index=index_seq)
 649     return seq
 650
 651 def get_unaligned_sample_export(lane, index_seq):
 652     body = """HWI-ST0787\t102\t{lane}\t1101\t1207\t1993\t{index}\t1\tAANGGATTCGATCCGGCTTAAGAGATGAAAACCGAAAGGGCCGACCGAA\taaBS`ccceg[`ae[dRR_[[SPPPP__ececfYYWaegh^\\ZLLY\\X`\tNM\t\t\t\t\t\t
 653 HWI-ST0787\t102     {lane}       1101    1478    1997    {index}  1       CAAGAACCCCGGGGGGGGGGGGGCAGAGAGGGGGAATTTTTTTTTTGTT       BBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBB       NM                                                                                      N
 654 HWI-ST0787      102     {lane}       1101    1625    1994    {index}  1       AANAATGCTACAGAGACAAAACAAAACTGATATGAAAGTTGAGAATAAA       \^BS\cccgegg[Q[QQQ[`egdgffbeggfgh^^YcfgfhXaHY^O^c       chrII.fa
 655 """.format(lane=lane, index=index_seq)
 656     return body
 657
 658 def ls_tree(root):
 659     for dirpath, dirnames, filenames in os.walk(root):
 660         for filename in filenames:
 661             print os.path.join(dirpath, filename)
 662
 663
 664 class BaseCallInfo(object):
 665     """Provide customization for how to setup the base call mock data
 666     """
 667     def __init__(self, qseq_file, tile_list, basecall_summary):
 668         self.qseq_file = qseq_file
 669         self.tile_list = tile_list
 670         self.basecall_summary = basecall_summary
 671
 672 # First generation HiSeq Flowcell
 673 ABXX_BASE_CALL_INFO = BaseCallInfo(
 674     qseq_file='AA01CCABXX_8_2_2207_qseq.txt',
 675     tile_list = HISEQ_TILE_LIST,
 676     basecall_summary = 'AA01CCABXX_BustardSummary.xml')