1 /* @source ajseqabi ***********************************************************
2 **
3 ** AJAX ABI format sequence processing functions
4 **
5 ** These functions control all aspects of AJAX ABI file processing
6 **
7 ** @author Copyright (C) 2000 Peter Rice
8 ** @version $Revision: 1.26 $
9 ** @modified 2000-2011 Peter Rice
10 ** @modified $Date: 2011/10/18 14:23:40 $ by $Author: rice $
11 ** @@
12 **
13 ** This library is free software; you can redistribute it and/or
14 ** modify it under the terms of the GNU Lesser General Public
15 ** License as published by the Free Software Foundation; either
16 ** version 2.1 of the License, or (at your option) any later version.
17 **
18 ** This library is distributed in the hope that it will be useful,
19 ** but WITHOUT ANY WARRANTY; without even the implied warranty of
20 ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
21 ** Lesser General Public License for more details.
22 **
23 ** You should have received a copy of the GNU Lesser General Public
24 ** License along with this library; if not, write to the Free Software
25 ** Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
26 ** MA 02110-1301, USA.
27 **
28 ******************************************************************************/
29
30 #include "ajlib.h"
31
32 #include "ajseqabi.h"
33 #include "ajfileio.h"
34
35
36
37
38 static AjBool seqABIReadInt4(AjPFile fp,ajlong *i4);
39 static AjBool seqABIReadFloat4(AjPFile fp,float* f4);
40 static AjBool seqABIReadInt2(AjPFile fp, ajshort *i2);
41 static AjBool seqABIGetFlag(AjPFile fp, ajlong flagLabel,
42 ajlong flagInstance, ajlong word, ajlong* val);
43 static AjBool seqABIGetFlagF(AjPFile fp, ajlong flagLabel,
44 ajlong flagInstance, ajlong word,float* val);
45 static AjBool seqABIGetFlagW(AjPFile fp, ajlong flagLabel,
46 ajlong word, ajshort* val);
47 static ajshort seqABIBaseIdx(char B);
48
49
50
51
52
53 /* @func ajSeqABITest *********************************************************
54 **
55 ** Test file type is ABI format - look for 'ABIF' flag (which may be in one
56 ** of 2 places).
57 **
58 ** @param [u] fp [AjPFile] ABI format file
59 ** @return [AjBool] ajTrue on success
60 **
61 ** @release 1.8.0
62 ** @@
63 ******************************************************************************/
64
ajSeqABITest(AjPFile fp)65 AjBool ajSeqABITest(AjPFile fp)
66 {
67 char pabi[5];
68 pabi[4] = '\0';
69
70 ajDebug("ajSeqABITest file %F end: %B\n", fp, fp->End);
71
72 if (fp->End && ajFileIsStdin(fp))
73 {
74 ajDebug("EOF: ajSeqABITest already at end file %F\n", fp);
75
76 return ajFalse;
77 }
78
79 if(ajFileSeek(fp,0,SEEK_SET) >= 0)
80 if(ajReadbinBinary(fp,1,4,(void *)pabi))
81 {
82 ajDebug("ajSeqABITest was at '%s'\n", pabi);
83
84 if(ajCharPrefixC(pabi,"ABIF"))
85 return ajTrue;
86 }
87
88 if(ajFileSeek(fp,26,SEEK_SET) >= 0)
89 {
90 ajDebug("ajSeqABITest seek to pos 26\n");
91
92 if(ajReadbinBinary(fp,1,4,(void*)pabi))
93 {
94 ajDebug("ajSeqABITest seek to '%s'\n", pabi);
95
96 if(ajCharPrefixC(pabi,"ABIF"))
97 return ajTrue;
98 }
99 }
100
101 return ajFalse;
102 }
103
104
105
106
107 /* @func ajSeqABIReadConfid ***************************************************
108 **
109 ** Read in confidence values from an ABI trace file.
110 **
111 ** @param [u] fp [AjPFile] ABI format input file
112 ** @param [r] pconO [ajlong] PCON offset in an ABI file
113 ** @param [r] numBases [ajlong] number of bases
114 ** @param [w] Pqual [float *] array of confidence values
115 ** @return [AjBool] ajTrue on success
116 **
117 ** @release 6.3.0
118 ** @@
119 ******************************************************************************/
120
ajSeqABIReadConfid(AjPFile fp,ajlong pconO,ajlong numBases,float * Pqual)121 AjBool ajSeqABIReadConfid(AjPFile fp,ajlong pconO,ajlong numBases,
122 float * Pqual)
123 {
124 ajint i;
125 char qv;
126
127 ajDebug("ajSeqABIReadConfid pcon0 %Ld numBases %Ld\n", pconO, numBases);
128
129 ajFileSeek(fp,pconO,SEEK_SET);
130
131 for (i=0;i<(ajint)numBases;i++)
132 {
133 ajReadbinBinary(fp,1,1,&qv);
134 Pqual[i] = (float)(int)qv;
135 }
136
137 return ajTrue;
138 }
139
140
141
142
143 /* @func ajSeqABIReadSeq ******************************************************
144 **
145 ** Read in a sequence from an ABI trace file.
146 **
147 ** @param [u] fp [AjPFile] ABI format input file
148 ** @param [r] baseO [ajlong] PBAS offset in an ABI file
149 ** @param [r] numBases [ajlong] number of bases
150 ** @param [w] nseq [AjPStr*] read sequence
151 ** @return [AjBool] ajTrue on success
152 **
153 ** @release 1.8.0
154 ** @@
155 ******************************************************************************/
156
ajSeqABIReadSeq(AjPFile fp,ajlong baseO,ajlong numBases,AjPStr * nseq)157 AjBool ajSeqABIReadSeq(AjPFile fp,ajlong baseO,ajlong numBases,
158 AjPStr* nseq)
159 {
160 ajint i;
161 char pseq;
162
163 ajDebug("ajSeqABIReadSeq base0 %Ld numBases %Ld\n", baseO, numBases);
164
165 ajFileSeek(fp,baseO,SEEK_SET);
166
167 for (i=0;i<(ajint)numBases;i++)
168 {
169 ajReadbinBinary(fp,1,1,&pseq);
170 ajStrAppendK(nseq,pseq);
171 }
172
173 return ajTrue;
174 }
175
176
177
178
179 /* @func ajSeqABIMachineName **************************************************
180 **
181 ** Get the name of the machine used to obtain an ABI trace file.
182 **
183 ** @param [u] fp [AjPFile] ABI format file
184 ** @param [w] machine [AjPStr*] machine name
185 ** @return [AjBool] ajTrue on success
186 **
187 ** @release 1.8.0
188 ** @@
189 ******************************************************************************/
190
ajSeqABIMachineName(AjPFile fp,AjPStr * machine)191 AjBool ajSeqABIMachineName(AjPFile fp,AjPStr *machine)
192 {
193 ajlong mchn;
194 ajlong MCHNtag;
195 unsigned char l;
196
197 MCHNtag = ((ajlong) ((((('M'<<8)+'C')<<8)+'H')<<8)+'N');
198 ajDebug("getflag MCHN\n");
199
200 if(seqABIGetFlag(fp,MCHNtag,1,5,&mchn))
201 {
202 if (ajFileSeek(fp,mchn,SEEK_SET) >= 0)
203 {
204 ajReadbinBinary(fp,1,sizeof(char),&l);
205 *machine = ajStrNewRes(l+1);
206 ajReadbinBinary(fp,1,l,(void*)ajStrGetuniquePtr(machine));
207 *(ajStrGetuniquePtr(machine)+l)='\0';
208 ajStrSetValid(machine);
209 }
210 else
211 return ajFalse;
212 }
213 else
214 return ajFalse;
215
216 return ajTrue;
217 }
218
219
220
221
222 /* @func ajSeqABIGetNData *****************************************************
223 **
224 ** Find 'DATA' tag and get the number of data points.
225 **
226 ** @param [u] fp [AjPFile] ABI format file
227 ** @return [ajint] Number of data points in file
228 **
229 ** @release 1.8.0
230 ** @@
231 ******************************************************************************/
232
ajSeqABIGetNData(AjPFile fp)233 ajint ajSeqABIGetNData(AjPFile fp)
234 {
235
236 ajlong numPoints;
237 ajlong DATAtag;
238 ajshort TRACE_INDEX;
239
240 DATAtag = ((ajlong) ((((('D'<<8)+'A')<<8)+'T')<<8)+'A');
241 TRACE_INDEX = 9;
242
243 if(!seqABIGetFlag(fp,DATAtag,TRACE_INDEX,3,&numPoints))
244 ajFatal("Error - locating DATA tag");
245
246 return (ajint) numPoints;
247 }
248
249
250
251
252 /* @func ajSeqABIGetNBase *****************************************************
253 **
254 ** Find the 'PBAS' tag in an ABI trace file and get the number of bases.
255 **
256 ** @param [u] fp [AjPFile] ABI format file
257 ** @return [ajint] Number of bases in file
258 **
259 ** @release 1.8.0
260 ** @@
261 ******************************************************************************/
262
ajSeqABIGetNBase(AjPFile fp)263 ajint ajSeqABIGetNBase(AjPFile fp)
264 {
265 ajlong numBases;
266 ajlong BASEtag;
267
268 BASEtag = ((ajlong) ((((('P'<<8)+'B')<<8)+'A')<<8)+'S');
269
270 ajDebug("getflag PBAS 2 NBase\n");
271 if(!seqABIGetFlag(fp,BASEtag,2,3,&numBases))
272 if(!seqABIGetFlag(fp,BASEtag,1,3,&numBases))
273 {
274 ajWarn("ABI file has no PBAS tag for sequence data");
275 return 0;
276 }
277
278 return (ajint) numBases;
279 }
280
281
282
283
284 /* @func ajSeqABIGetData ******************************************************
285 **
286 ** Read in the processed trace data from an ABI file.
287 **
288 ** @param [u] fp [AjPFile] ABI format file
289 ** @param [r] Offset [const ajlong*] data offset in ABI file
290 ** @param [r] numPoints [ajlong] number of data points
291 ** @param [w] trace [AjPInt2d] (4xnumPoints) array of trace data
292 ** @return [void]
293 **
294 ** @release 1.8.0
295 ** @@
296 ******************************************************************************/
297
ajSeqABIGetData(AjPFile fp,const ajlong * Offset,ajlong numPoints,AjPInt2d trace)298 void ajSeqABIGetData(AjPFile fp,const ajlong *Offset,ajlong numPoints,
299 AjPInt2d trace)
300 {
301 ajint i;
302 ajint j;
303 ajshort traceValue;
304
305 /* Read in data */
306 for (i=0;i<4;i++)
307 {
308 if (ajFileSeek(fp,Offset[i],SEEK_SET))
309 ajFatal("Error - reading trace");
310
311 for (j=0;j<(ajint)numPoints;j++)
312 if (seqABIReadInt2(fp,&traceValue))
313 ajInt2dPut(&trace,i,j,(ajint)traceValue);
314 else
315 ajFatal("Error - reading trace");
316 }
317
318 return;
319 }
320
321
322
323
324 /* @func ajSeqABIGetBasePosition **********************************************
325 **
326 ** Read in the base positions from an ABI file.
327 **
328 ** @param [u] fp [AjPFile] ABI format file
329 ** @param [r] numBases [ajlong] number of bases to be read
330 ** @param [w] basePositions [AjPShort*] base positions output
331 ** @return [void]
332 **
333 ** @release 1.8.0
334 ** @@
335 ******************************************************************************/
336
ajSeqABIGetBasePosition(AjPFile fp,ajlong numBases,AjPShort * basePositions)337 void ajSeqABIGetBasePosition(AjPFile fp,ajlong numBases,
338 AjPShort* basePositions)
339 {
340 ajint i;
341 ajshort bP;
342
343 /* Read in base positions */
344 for (i=0;i<(ajint)numBases;i++)
345 {
346 if (!seqABIReadInt2(fp,&bP))
347 ajFatal("Error - in finding Base Position");
348
349 ajShortPut(basePositions,i,bP);
350 }
351
352 return;
353 }
354
355
356
357
358 /* @func ajSeqABIGetSignal ****************************************************
359 **
360 ** Read in the signal strength information from an ABI file.
361 **
362 ** @param [u] fp [AjPFile] ABI format file
363 ** @param [r] fwo_ [ajlong] field order
364 ** @param [w] sigC [ajshort*] average signal strength for C
365 ** @param [w] sigA [ajshort*] average signal strength for A
366 ** @param [w] sigG [ajshort*] average signal strength for G
367 ** @param [w] sigT [ajshort*] average signal strength for T
368 ** @return [void]
369 **
370 ** @release 1.8.0
371 ** @@
372 ******************************************************************************/
373
ajSeqABIGetSignal(AjPFile fp,ajlong fwo_,ajshort * sigC,ajshort * sigA,ajshort * sigG,ajshort * sigT)374 void ajSeqABIGetSignal(AjPFile fp,ajlong fwo_,
375 ajshort *sigC,ajshort *sigA,
376 ajshort *sigG,ajshort *sigT)
377 {
378 ajlong signalO;
379 ajshort* base[4];
380
381 ajlong SIGNALtag;
382
383 SIGNALtag = ((ajlong) ((((('S'<<8)+'/')<<8)+'N')<<8)+'%');
384 ajDebug("getflag S/N%%\n");
385
386 /* Get signal strength info */
387 if (seqABIGetFlag(fp,SIGNALtag,1,5,&signalO))
388 {
389 base[0] = sigC;
390 base[1] = sigA;
391 base[2] = sigG;
392 base[3] = sigT;
393 if (ajFileSeek(fp, signalO, SEEK_SET) >= 0 &&
394 seqABIReadInt2(fp, base[seqABIBaseIdx((char)(fwo_>>24&255))]) &&
395 seqABIReadInt2(fp, base[seqABIBaseIdx((char)(fwo_>>16&255))]) &&
396 seqABIReadInt2(fp, base[seqABIBaseIdx((char)(fwo_>>8&255))]) &&
397 seqABIReadInt2(fp, base[seqABIBaseIdx((char)(fwo_&255))]))
398 {
399 /*
400 ajDebug("avg_signal_strength = C:%d A:%d G:%d T:%d\n",sigC,sigA,
401 sigG,sigT);
402 */
403 }
404 }
405
406 return;
407 }
408
409
410
411
412 /* @func ajSeqABIGetBaseSpace *************************************************
413 **
414 ** Read in the base spacing from an ABI file.
415 **
416 ** @param [u] fp [AjPFile] ABI format file
417 ** @return [float] base spacing
418 **
419 ** @release 1.8.0
420 ** @@
421 ******************************************************************************/
422
ajSeqABIGetBaseSpace(AjPFile fp)423 float ajSeqABIGetBaseSpace(AjPFile fp)
424 {
425
426 float spacing = 0.;
427 ajlong SPACINGtag;
428
429 SPACINGtag = ((ajlong) ((((('S'<<8)+'P')<<8)+'A')<<8)+'C');
430
431 seqABIGetFlagF(fp,SPACINGtag,1,5,&spacing);
432
433 return spacing;
434 }
435
436
437
438
439 /* @func ajSeqABIGetBaseOffset ************************************************
440 **
441 ** Routine to get the 'PBAS' tag offset in an ABI file.
442 **
443 ** @param [u] fp [AjPFile] ABI format file
444 ** @return [ajint] 'PBAS' tag offset in an ABI file
445 **
446 ** @release 1.8.0
447 ** @@
448 ******************************************************************************/
449
ajSeqABIGetBaseOffset(AjPFile fp)450 ajint ajSeqABIGetBaseOffset(AjPFile fp)
451 {
452 ajlong baseO;
453 ajlong BASEtag;
454
455 BASEtag = ((ajlong) ((((('P'<<8)+'B')<<8)+'A')<<8)+'S');
456 ajDebug("getflag PBAS 2 offset\n");
457
458 /* Find BASE tag & get offset */
459 if(!seqABIGetFlag(fp,BASEtag,2,5,&baseO))
460 if(!seqABIGetFlag(fp,BASEtag,1,5,&baseO))
461 {
462 ajWarn("ABI file has no PBAS flag for Base Offset");
463 return 0;
464 }
465
466 return (ajint) baseO;
467 }
468
469
470
471
472 /* @func ajSeqABIGetBasePosOffset *********************************************
473 **
474 ** Routine to get the 'PLOC', base position, tag offset in an ABI file
475 **
476 ** @param [u] fp [AjPFile] ABI format file
477 ** @return [ajint] base position offset in an ABI file
478 **
479 ** @release 1.8.0
480 ** @@
481 ******************************************************************************/
482
ajSeqABIGetBasePosOffset(AjPFile fp)483 ajint ajSeqABIGetBasePosOffset(AjPFile fp)
484 {
485 ajlong basePosO;
486 ajlong BASEPOStag;
487
488 BASEPOStag = ((ajlong) ((((('P'<<8)+'L')<<8)+'O')<<8)+'C');
489 ajDebug("getflag PLOC 2\n");
490
491 /* Find BASEPOS tag & get base position offset */
492 if (!seqABIGetFlag(fp,BASEPOStag,2,5,&basePosO))
493 if (!seqABIGetFlag(fp,BASEPOStag,1,5,&basePosO))
494 {
495 ajWarn("ABI file has no PLOC flag for Base Pos Offset");
496 return 0;
497 }
498
499 return (ajint) basePosO;
500 }
501
502
503
504
505 /* @func ajSeqABIGetConfidOffset **********************************************
506 **
507 ** Routine to get the 'PCON' confidence tag offset in an ABI file.
508 **
509 ** @param [u] fp [AjPFile] ABI format file
510 ** @return [ajint] 'PCON' tag offset in an ABI file
511 **
512 ** @release 6.3.0
513 ** @@
514 ******************************************************************************/
515
ajSeqABIGetConfidOffset(AjPFile fp)516 ajint ajSeqABIGetConfidOffset(AjPFile fp)
517 {
518 ajlong pconO;
519 ajlong PCONtag;
520
521 PCONtag = ((ajlong) ((((('P'<<8)+'C')<<8)+'O')<<8)+'N');
522 ajDebug("getflag PCON 2\n");
523
524 /* Find PCON tag & get offset */
525 if(!seqABIGetFlag(fp,PCONtag,2,5,&pconO))
526 if(!seqABIGetFlag(fp,PCONtag,1,5,&pconO))
527 {
528 /* ajWarn("ABI file has no PCON flag for Confidence Offset");*/
529 return 0;
530 }
531
532 return (ajint) pconO;
533 }
534
535
536
537
538 /* @func ajSeqABIGetFWO *******************************************************
539 **
540 ** Routine to get the "FWO" tag, field order ("GATC"), tag.
541 **
542 ** @param [u] fp [AjPFile] ABI format file
543 ** @return [ajint] field order
544 **
545 ** @release 1.8.0
546 ** @@
547 ******************************************************************************/
548
ajSeqABIGetFWO(AjPFile fp)549 ajint ajSeqABIGetFWO(AjPFile fp)
550 {
551
552 ajlong fwo_;
553 ajlong FWO_tag;
554
555 FWO_tag = ((ajlong) ((((('F'<<8)+'W')<<8)+'O')<<8)+'_');
556 ajDebug("getflag FWO_ 2\n");
557
558 /* Find FWO tag */
559 if (!seqABIGetFlag(fp,FWO_tag,2,5,&fwo_))
560 if (!seqABIGetFlag(fp,FWO_tag,1,5,&fwo_))
561 {
562 ajWarn("ABI file has no FWO_ flag for field order");
563 return 0;
564 }
565
566 return (ajint) fwo_;
567 }
568
569
570
571
572 /* @func ajSeqABIGetPrimerOffset **********************************************
573 **
574 ** Routine to get the primer offset in an ABI file.
575 **
576 ** @param [u] fp [AjPFile] ABI format file
577 ** @return [ajint] primer offset
578 **
579 ** @release 1.8.0
580 ** @@
581 ******************************************************************************/
582
ajSeqABIGetPrimerOffset(AjPFile fp)583 ajint ajSeqABIGetPrimerOffset(AjPFile fp)
584 {
585 ajshort primerPos;
586 ajlong PPOStag;
587
588 PPOStag = ((ajlong) ((((('P'<<8)+'P')<<8)+'O')<<8)+'S');
589
590 /* Find PPOS tag (Primer Position) & get offset */
591 if (!seqABIGetFlagW(fp,PPOStag,6,&primerPos))
592 {
593 ajWarn("ABI file has no PPOS flag for primer offset");
594 return 0;
595 }
596
597 return primerPos;
598 }
599
600
601
602
603 /* @func ajSeqABIGetPrimerPosition ********************************************
604 **
605 ** Routine to get the primer position in an ABI file.
606 **
607 ** @param [u] fp [AjPFile] ABI format file
608 ** @return [ajint] primer position
609 **
610 ** @release 1.8.0
611 ** @@
612 ******************************************************************************/
613
ajSeqABIGetPrimerPosition(AjPFile fp)614 ajint ajSeqABIGetPrimerPosition(AjPFile fp)
615 {
616 ajlong primerPosition;
617 ajlong PPOStag;
618
619 PPOStag = ((ajlong) ((((('P'<<8)+'P')<<8)+'O')<<8)+'S');
620 ajDebug("getflag PPOS 2\n");
621
622 if (!seqABIGetFlag(fp,PPOStag,2,5,&primerPosition))
623 if (!seqABIGetFlag(fp,PPOStag,1,5,&primerPosition))
624 {
625 ajWarn("ABI file has no PPOS flag for primer position");
626 return 0;
627 }
628
629 /* ppos stored in MBShort of pointer */
630 primerPosition = primerPosition>>16;
631
632
633 return (ajint) primerPosition;
634 }
635
636
637
638
639 /* @func ajSeqABIGetTraceOffset ***********************************************
640 **
641 ** Get the processed trace data ('DATA' tag) offset in an ABI file.
642 **
643 ** @param [u] fp [AjPFile] ABI format file
644 ** @param [w] Offset [ajlong *] trace data offset, used in ajSeqABIGetData
645 ** @return [AjBool] ajTrue on success
646 **
647 ** @release 1.8.0
648 ** @@
649 ******************************************************************************/
650
ajSeqABIGetTraceOffset(AjPFile fp,ajlong * Offset)651 AjBool ajSeqABIGetTraceOffset(AjPFile fp, ajlong *Offset)
652 {
653 ajlong dataxO[4];
654 ajlong fwo_;
655
656 /* BYTE[i] is a byte mask for byte i */
657 const ajlong BYTE[] = { 0x000000ff };
658 ajshort TRACE_INDEX;
659 ajlong DATAtag;
660
661 TRACE_INDEX = 9;
662 DATAtag = ((ajlong) ((((('D'<<8)+'A')<<8)+'T')<<8)+'A');
663
664 /* Find FWO tag - Field order "GATC" */
665 fwo_ = ajSeqABIGetFWO(fp);
666
667 /* Get data trace offsets */
668 if (!seqABIGetFlag(fp,DATAtag,TRACE_INDEX,
669 5,&dataxO[seqABIBaseIdx((char)(fwo_>>24&BYTE[0]))]))
670 return ajFalse;
671
672 if (!seqABIGetFlag(fp,DATAtag,TRACE_INDEX+1,
673 5,&dataxO[seqABIBaseIdx((char)(fwo_>>16&BYTE[0]))]))
674 return ajFalse;
675
676 if (!seqABIGetFlag(fp,DATAtag,TRACE_INDEX+2,
677 5,&dataxO[seqABIBaseIdx((char)(fwo_>>8&BYTE[0]))]))
678 return ajFalse;
679
680 if (!seqABIGetFlag(fp,DATAtag,TRACE_INDEX+3,
681 5,&dataxO[seqABIBaseIdx((char)(fwo_&BYTE[0]))]))
682 return ajFalse;
683
684 Offset[0]=dataxO[seqABIBaseIdx((char)(fwo_>>24&BYTE[0]))];
685 Offset[1]=dataxO[seqABIBaseIdx((char)(fwo_>>16&BYTE[0]))];
686 Offset[2]=dataxO[seqABIBaseIdx((char)(fwo_>>8&BYTE[0]))];
687 Offset[3]=dataxO[seqABIBaseIdx((char)(fwo_&BYTE[0]))];
688
689 return ajTrue;
690 }
691
692
693
694
695 /* @funcstatic seqABIReadInt4 ************************************************
696 **
697 ** Routine to read 4 bytes from a file and return the integer.
698 **
699 ** @param [u] fp [AjPFile] ABI format file
700 ** @param [w] i4 [ajlong *] ajlong integer read in from ABI file
701 ** @return [AjBool] true if read successfully
702 **
703 ** @release 2.0.0
704 ** @@
705 ******************************************************************************/
706
seqABIReadInt4(AjPFile fp,ajlong * i4)707 static AjBool seqABIReadInt4(AjPFile fp,ajlong *i4)
708 {
709
710 unsigned char buf[sizeof(ajlong)];
711
712 if (ajReadbinBinary(fp,1,4,(void *)buf) != 1)
713 return ajFalse;
714
715 *i4 = (ajlong)
716 (((ajulong)buf[3]) +
717 ((ajulong)buf[2]<<8) +
718 ((ajulong)buf[1]<<16) +
719 ((ajulong)buf[0]<<24));
720
721 /*ajDebug("seqABIReadInt4 %Ld, %x %x %x %x\n", *i4,
722 (int) buf[0],(int) buf[1],
723 (int) buf[2],(int) buf[3]);*/
724
725 return ajTrue;
726 }
727
728
729
730
731 /* @funcstatic seqABIReadFloat4 ***********************************************
732 **
733 ** Routine to read 4 bytes from a file and return the float.
734 **
735 ** @param [u] fp [AjPFile] ABI format file
736 ** @param [w] f4 [float *] float read in from ABI file
737 ** @return [AjBool] true if read successfully
738 **
739 ** @release 1.8.0
740 ** @@
741 ******************************************************************************/
742
seqABIReadFloat4(AjPFile fp,float * f4)743 static AjBool seqABIReadFloat4(AjPFile fp,float* f4)
744 {
745
746 unsigned char buf[sizeof(ajlong)];
747 ajulong res;
748
749 if (ajReadbinBinary(fp,1,4,(void *)buf) != 1)
750 return ajFalse;
751
752 res = (ajulong)
753 (((ajulong)buf[3]) +
754 ((ajulong)buf[2]<<8) +
755 ((ajulong)buf[1]<<16) +
756 ((ajulong)buf[0]<<24));
757
758 *f4 = (float) res;
759
760 return ajTrue;
761 }
762
763
764
765
766 /* @funcstatic seqABIReadInt2 *************************************************
767 **
768 ** Routine to read 2 bytes from a file and return the short integer.
769 **
770 ** @param [u] fp [AjPFile] ABI format file
771 ** @param [w] i2 [ajshort *] short integer read in from ABI file
772 ** @return [AjBool] true if read successfully
773 **
774 ** @release 1.8.0
775 ** @@
776 ******************************************************************************/
777
seqABIReadInt2(AjPFile fp,ajshort * i2)778 static AjBool seqABIReadInt2(AjPFile fp, ajshort *i2)
779 {
780 unsigned char buf[sizeof(ajshort)];
781
782 if (ajReadbinBinary(fp,1,2,(void *)buf) != 1)
783 return ajFalse;
784
785 *i2 = (ajshort)
786 (((ajushort)buf[1]) +
787 ((ajushort)buf[0]<<8));
788
789 return ajTrue;
790 }
791
792
793
794
795 /* @funcstatic seqABIGetFlag **************************************************
796 **
797 ** Routine to read through an ABI trace file until it reaches a flag
798 ** (flagLabel). If there are multiple flags in the file it will search
799 ** to find the correct instance of that flag (flagInstance).
800 ** It will then return the *integer* value (val) of the word+1 from
801 ** that flag record.
802 **
803 ** @param [u] fp [AjPFile] ABI format file
804 ** @param [r] flagLabel [ajlong] flag in the ABI file
805 ** @param [r] flagInstance [ajlong] flag instance in the ABI file
806 ** @param [r] word [ajlong] number of fields to ignore in this record
807 ** @param [w] val [ajlong*] integer value of the word+1
808 ** @return [AjBool] true if read successfully
809 **
810 ** @release 1.8.0
811 ** @@
812 ******************************************************************************/
813
seqABIGetFlag(AjPFile fp,ajlong flagLabel,ajlong flagInstance,ajlong word,ajlong * val)814 static AjBool seqABIGetFlag(AjPFile fp, ajlong flagLabel,
815 ajlong flagInstance, ajlong word, ajlong* val)
816 {
817 ajint flagNum = -1;
818 ajint i;
819 ajlong Label;
820 ajlong Instance;
821 ajlong indexO;
822 ajint INDEX_ENTRY_LENGTH;
823
824 INDEX_ENTRY_LENGTH= 28;
825
826 if(ajFileSeek(fp,26,SEEK_SET) ||
827 (!seqABIReadInt4(fp, &indexO))) ajFatal("Error - in finding flag");
828
829 ajDebug("seqABIGetFlag Flag %Lp %Ld %Ld indexO %Ld\n",
830 flagLabel, flagInstance, word, indexO);
831 do
832 {
833 flagNum++;
834
835 if (ajFileSeek(fp,indexO+(flagNum*INDEX_ENTRY_LENGTH),SEEK_SET) != 0)
836 {
837 ajDebug("num %d Flag %Lp %Ld %Ld end of file\n",
838 flagNum, flagLabel, flagInstance, word);
839 return ajFalse;
840 }
841
842 if (!seqABIReadInt4(fp, &Label))
843 return ajFalse;
844
845 if (!seqABIReadInt4(fp, &Instance))
846 return ajFalse;
847
848 ajDebug("Label %Lp Instance %Ld\n", Label, Instance);
849
850 } while (!(Label == (ajlong)flagLabel &&
851 Instance == (ajlong)flagInstance));
852
853 for (i=2; i<=word; i++)
854 if (!seqABIReadInt4(fp, val))
855 return ajFalse;
856
857 ajDebug("OK Flag %Lp %Ld %Ld (%d)\n",
858 flagLabel, flagInstance, word, val[0]);
859 return ajTrue;
860
861 }
862
863
864
865
866 /* @funcstatic seqABIGetFlagF *************************************************
867 **
868 ** Routine to read through an ABI trace file until it reaches a flag
869 ** (flagLabel). If there are multiple flags in the file it will search
870 ** to find the correct instance of that flag (flagInstance).
871 ** It will then return the *float* value (val) of the word+1 from
872 ** that flag record.
873 **
874 ** @param [u] fp [AjPFile] ABI format file
875 ** @param [r] flagLabel [ajlong] flag in the ABI file
876 ** @param [r] flagInstance [ajlong] flag instance in the ABI file
877 ** @param [r] word [ajlong] number of fields to ignore in this record
878 ** @param [w] val [float*] integer value of the word+1
879 ** @return [AjBool] true if read successfully
880 **
881 ** @release 1.8.0
882 ** @@
883 ******************************************************************************/
884
seqABIGetFlagF(AjPFile fp,ajlong flagLabel,ajlong flagInstance,ajlong word,float * val)885 static AjBool seqABIGetFlagF(AjPFile fp, ajlong flagLabel,
886 ajlong flagInstance, ajlong word,float* val)
887 {
888 ajint flagNum = -1;
889 ajint i;
890 ajlong Label;
891 ajlong Instance;
892 ajlong indexO;
893 ajint INDEX_ENTRY_LENGTH;
894
895 INDEX_ENTRY_LENGTH= 28;
896
897 if(ajFileSeek(fp,26,SEEK_SET) ||
898 (!seqABIReadInt4(fp, &indexO))) ajFatal("Error - in finding flag");
899
900 do
901 {
902 flagNum++;
903
904 if (ajFileSeek(fp,indexO+(flagNum*INDEX_ENTRY_LENGTH),SEEK_SET) != 0)
905 return ajFalse;
906
907 if (!seqABIReadInt4(fp, &Label))
908 return ajFalse;
909
910 if (!seqABIReadInt4(fp, &Instance))
911 return ajFalse;
912 } while (!(Label == (ajlong)flagLabel &&
913 Instance == (ajlong)flagInstance));
914
915 for (i=2; i<=word; i++)
916 if (!seqABIReadFloat4(fp, val))
917 return ajFalse;
918
919 return ajTrue;
920 }
921
922
923
924
925 /* @funcstatic seqABIGetFlagW *************************************************
926 **
927 ** Routine to read through an ABI trace file until it reaches a flag
928 ** (flagLabel). If there are multiple flags in the file it will search
929 ** to find the correct instance of that flag (flagInstance).
930 ** It will then return the *short ajint* value (val) of the word+1 from
931 ** that flag record.
932 **
933 ** @param [u] fp [AjPFile] ABI format file
934 ** @param [r] flagLabel [ajlong] flag in the ABI file
935 ** @param [r] word [ajlong] number of fields to ignore in this record
936 ** @param [w] val [ajshort*] integer value of the word+1
937 ** @return [AjBool] true if read successfully
938 **
939 ** @release 1.13.0
940 ** @@
941 ******************************************************************************/
942
seqABIGetFlagW(AjPFile fp,ajlong flagLabel,ajlong word,ajshort * val)943 static AjBool seqABIGetFlagW(AjPFile fp, ajlong flagLabel,
944 ajlong word, ajshort* val)
945 {
946 ajint flagNum = -1;
947 ajint i;
948 ajlong Label;
949 ajlong jval;
950 ajlong indexO;
951 ajint INDEX_ENTRY_LENGTH;
952
953 INDEX_ENTRY_LENGTH= 28;
954
955 if(ajFileSeek(fp,26,SEEK_SET) ||
956 (!seqABIReadInt4(fp, &indexO))) ajFatal("Error - in finding flag");
957
958 do
959 {
960 flagNum++;
961
962 if (ajFileSeek(fp, indexO+(flagNum*INDEX_ENTRY_LENGTH), SEEK_SET) != 0)
963 return ajFalse;
964
965 if (!seqABIReadInt4(fp, &Label))
966 return ajFalse;
967 }
968 while (Label != (ajlong)flagLabel);
969
970
971 for (i=2; i<word; i++)
972 if (!seqABIReadInt4(fp, &jval))
973 return ajFalse;
974
975 if (!seqABIReadInt2(fp, val))
976 return ajFalse;
977
978 return ajTrue;
979 }
980
981
982
983
984 /* @funcstatic seqABIBaseIdx **************************************************
985 **
986 ** Returns: 0 if C, 1 if A, 2 if G, 3 if anything else
987 **
988 ** @param [r] B [char] base (C, A, G or T)
989 ** @return [ajshort] 0 if C, 1 if A, 2 if G, 3 if anything else
990 **
991 ** @release 1.8.0
992 ** @@
993 ******************************************************************************/
994
seqABIBaseIdx(char B)995 static ajshort seqABIBaseIdx(char B)
996 {
997 return ((B)=='C' ? 0 : (B)=='A' ? 1 : (B)=='G' ? 2 : 3);
998 }
999
1000
1001
1002
1003 /* @func ajSeqABISampleName ***************************************************
1004 **
1005 ** Get the sample name from an ABI trace file.
1006 **
1007 ** @param [u] fp [AjPFile] ABI format file
1008 ** @param [w] sample [AjPStr*] sample name
1009 ** @return [AjBool] true if read successfully
1010 **
1011 ** @release 2.0.0
1012 ** @@
1013 ******************************************************************************/
1014
ajSeqABISampleName(AjPFile fp,AjPStr * sample)1015 AjBool ajSeqABISampleName(AjPFile fp, AjPStr *sample)
1016 {
1017 ajlong mchn;
1018 ajlong SMPLtag;
1019 unsigned char l;
1020
1021 SMPLtag = ((ajlong) ((((('S'<<8)+'M')<<8)+'P')<<8)+'L');
1022 ajDebug("getflag SMPL\n");
1023
1024 if((seqABIGetFlag(fp,SMPLtag,1,5,&mchn)) &&
1025 (ajFileSeek(fp,mchn,SEEK_SET) >= 0))
1026 {
1027 ajReadbinBinary(fp,1,sizeof(char),&l);
1028 *sample = ajStrNewRes(l+1);
1029 ajReadbinBinary(fp,1,l,(void*)ajStrGetuniquePtr(sample));
1030 *(ajStrGetuniquePtr(sample)+l)='\0';
1031 ajStrSetValid(sample);
1032 ajDebug("read SMPL at %Ld %u '%S'\n", mchn, (ajuint) l, *sample);
1033 }
1034
1035 return ajTrue;
1036 }
1037