1 /*****************************************************************************
2 *
3 * seqtest.c
4 * test program for sequence display, SeqPort and ToFasta
5 *
6 *****************************************************************************/
7
8 #include <seqport.h>
9 #include <tofasta.h>
10
11 void BuildList (SeqEntryPtr sep, Pointer data, Int4 index, Int2 indent);
12
13 #define NUMARG 2
14 Args myargs[NUMARG] = {
15 {"Filename for asn.1 input","stdin",NULL,NULL,FALSE,'i',ARG_FILE_IN,0.0,0,NULL},
16 {"Filename for output","stdout", NULL,NULL,TRUE,'o',ARG_FILE_OUT,0.0,0,NULL}};
17
Main(void)18 Int2 Main(void)
19 {
20 AsnIoPtr aip;
21 SeqEntryPtr sep;
22 BioseqPtr PNTR seqlist;
23 Int4 seqnum, i, numseg, lens[10], j;
24 Int2 ctr;
25 SeqPortPtr spp;
26 Uint1 residue;
27 FILE* fp;
28 CharPtr title;
29 Char buffer[101];
30 MonitorPtr mon;
31
32 /* check command line arguments */
33
34 if ( ! GetArgs("SeqTest",NUMARG, myargs))
35 return 1;
36
37 mon = MonitorStrNew("SeqTest", 40);
38 SetProgMon(StdProgMon, (Pointer)mon);
39
40 /*
41 ** Load SeqEntry object loader and sequence alphabets
42 */
43
44 if (! SeqEntryLoad()) {
45 Message(MSG_ERROR, "SeqEntryLoad failed");
46 return 1;
47 }
48
49 /*
50 ** Use the file "example.prt" as the ASN I/O stream. This file
51 ** can be found in the ncbi/demo. It is in ASN.1 Print Value format.
52 */
53
54 if ((aip = AsnIoOpen(myargs[0].strvalue, "r")) == NULL)
55 return 1;
56
57 /*
58 ** Write the output to "seqtest.out".
59 */
60
61 fp = FileOpen(myargs[1].strvalue, "w");
62 fprintf(fp, "Sequence summary:\n\n");
63
64 /*
65 ** Read in the whole entry into the Sequence Entry Pointer, sep.
66 ** Close the ASN stream, which in turn closes the input file.
67 */
68
69 sep = SeqEntryAsnRead(aip, NULL);
70 aip = AsnIoClose(aip);
71
72 mon = MonitorFree(mon);
73 SetProgMon(NULL, NULL);
74
75 /*
76 ** Determine how many Bioseqs are in this SeqEntry. Allocate
77 ** enough memory to hold a list of pointer to all of these
78 ** Bioseqs. Invoke an Explore function to "visit"each Bioseq.
79 ** We are allowed to pass one pointer for use by the exploring
80 ** function, in this case, "BuildList".
81 */
82
83 seqnum = BioseqCount(sep);
84 seqlist = MemNew((size_t)(seqnum * sizeof(BioseqPtr)));
85 BioseqExplore(sep, (Pointer) seqlist, BuildList);
86
87 /*
88 ** For each Bioseq in the SeqEntry write out it's title
89 ** len, number of gaps, and number of segments. Write out
90 ** the length of each segment, up to 10.
91 */
92
93 for(i = 0; i < seqnum; i++) {
94 numseg = BioseqCountSegs(seqlist[i]);
95 title = BioseqGetTitle(seqlist[i]);
96 FilePuts((VoidPtr)title, fp);
97 FilePuts("\n", fp);
98 fprintf(fp, "len=%d gaps=%d segs=%d\n", BioseqGetLen(seqlist[i]),
99 BioseqGetGaps(seqlist[i]), numseg);
100 if ((numseg > 1) && (numseg <= 10)) {
101 BioseqGetSegLens (seqlist[i], lens);
102 for (j = 0; j < numseg; j++)
103 fprintf(fp, " len = %d\n", lens[j]);
104 }
105 FilePuts("\n", fp);
106 }
107
108 spp = SeqPortNew(seqlist[0], 0, -1, 0, Seq_code_iupacna);
109 if (spp == NULL)
110 Message(MSG_ERROR, "fail on SeqPortNew");
111
112 fprintf(fp, "SeqPort: plus strand with SeqPortGetResidue\n\n");
113
114 i = 0;
115 while ((residue = SeqPortGetResidue(spp)) != SEQPORT_EOF) {
116 if (! IS_residue(residue)) {
117 buffer[i] = '\0';
118 fprintf(fp, "%s\n", buffer);
119 i = 0;
120 switch (residue)
121 {
122 case SEQPORT_VIRT:
123 fprintf(fp, "[Gap]\n");
124 break;
125 case SEQPORT_EOS:
126 fprintf(fp, "[EOS]\n");
127 break;
128 default:
129 fprintf(fp, "[Invalid Residue]\n");
130 break;
131 }
132
133 }
134 else {
135 buffer[i] = residue;
136 i++;
137 if (i == 60) {
138 buffer[i] = '\0';
139 fprintf(fp, "%s\n", buffer);
140 i = 0;
141 }
142 }
143 }
144
145 if (i) {
146 buffer[i] = '\0';
147 fprintf(fp, "%s\n", buffer);
148 }
149
150 fprintf(fp, "[EOF]\n");
151 SeqPortFree(spp);
152
153 fprintf(fp, "\nSeqPort on minus with SeqPortRead\n\n");
154 spp = SeqPortNew(seqlist[0], 0, -1, Seq_strand_minus, Seq_code_iupacna);
155
156 if (spp == NULL)
157 Message(MSG_ERROR, "fail on SeqPortNew");
158
159 do {
160 ctr = SeqPortRead(spp, (Uint1Ptr)buffer, 60);
161 if (ctr > 0) {
162 buffer[ctr] = '\0';
163 fprintf(fp, "%s\n", buffer);
164 } else {
165 ctr *= -1;
166 switch (ctr)
167 {
168 case SEQPORT_VIRT:
169 fprintf(fp, "[Gap]\n");
170 break;
171 case SEQPORT_EOS:
172 fprintf(fp, "[EOS]\n");
173 break;
174 case SEQPORT_EOF:
175 fprintf(fp, "[EOF]\n");
176 break;
177 default:
178 fprintf(fp, "[Invalid Residue]\n");
179 break;
180 }
181 }
182 } while (ctr != SEQPORT_EOF);
183
184 SeqPortFree(spp);
185
186 /*
187 ** Write out the nucleic acid sequences in this SeqEntry
188 */
189
190 fprintf(fp, "\nNucleic Acids in FASTA format:\n\n");
191 SeqEntryToFasta(sep, fp, TRUE);
192
193 /*
194 ** Write out the protein sequences in this SeqEntry.
195 */
196
197 fprintf(fp, "\nProteins in FASTA format:\n\n");
198 SeqEntryToFasta(sep, fp, FALSE);
199
200 /*
201 ** Close the output file and free up allocated space.
202 */
203
204 fclose(fp);
205 MemFree(seqlist);
206 SeqEntryFree(sep);
207
208 return 0;
209 }
210
211 /*
212 ** This SeqEntry exploration function copy the current pointer position inthe
213 ** the Bioseq entry to a list of Bioseq pointers
214 */
215
BuildList(SeqEntryPtr sep,Pointer data,Int4 index,Int2 indent)216 void BuildList(SeqEntryPtr sep, Pointer data, Int4 index, Int2 indent)
217 {
218 ((BioseqPtr PNTR) data)[index] = (BioseqPtr)sep->data.ptrvalue;
219 return;
220 }
221
222
223
224