1 /* @source dottup application
2 **
3 ** Dotplot of two sequences
4 **
5 ** @author Copyright (C) Ian Longden
6 ** @modified: Alan Bleasby. Added non-proportional plot
7 ** @@
8 **
9 ** This program is free software; you can redistribute it and/or
10 ** modify it under the terms of the GNU General Public License
11 ** as published by the Free Software Foundation; either version 2
12 ** of the License, or (at your option) any later version.
13 **
14 ** This program is distributed in the hope that it will be useful,
15 ** but WITHOUT ANY WARRANTY; without even the implied warranty of
16 ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 ** GNU General Public License for more details.
18 **
19 ** You should have received a copy of the GNU General Public License
20 ** along with this program; if not, write to the Free Software
21 ** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
22 ******************************************************************************/
23
24 #include "emboss.h"
25
26
27
28
29 static void dottup_drawPlotlines(void *x, void *cl);
30 static void dottup_plotMatches(const AjPList list);
31 static void dottup_stretchplot(AjPGraph graph, const AjPList matchlist,
32 const AjPSeq seq1, const AjPSeq seq2,
33 ajint begin1, ajint begin2,
34 ajint end1, ajint end2);
35
36 static PLFLT offset1;
37 static PLFLT offset2;
38
39
40
41
42 /* @prog dottup ***************************************************************
43 **
44 ** Displays a wordmatch dotplot of two sequences
45 **
46 ******************************************************************************/
47
main(int argc,char ** argv)48 int main(int argc, char **argv)
49 {
50 AjPSeq seq1;
51 AjPSeq seq2;
52 ajint wordlen;
53 AjPTable seq1MatchTable = 0;
54 AjPList matchlist = NULL;
55 AjPGraph graph = NULL;
56 AjPGraph xygraph = NULL;
57 AjBool boxit;
58 /*
59 ** Different ticks as they need to be different for x and y due to
60 ** length of string being important on x
61 */
62 ajuint acceptableticksx[]=
63 {
64 1,10,50,100,500,1000,1500,10000,
65 500000,1000000,5000000
66 };
67 ajuint acceptableticks[]=
68 {
69 1,10,50,100,200,500,1000,2000,5000,10000,15000,
70 500000,1000000,5000000
71 };
72 ajint numbofticks = 10;
73 float xmargin;
74 float ymargin;
75 float ticklen;
76 float tickgap;
77 float onefifth = 0.0;
78 ajint i;
79 float k;
80 float max;
81 char ptr[10];
82 ajint begin1;
83 ajint begin2;
84 ajint end1;
85 ajint end2;
86 ajuint len1;
87 ajuint len2;
88 float fbegin1;
89 float fbegin2;
90 float fend1;
91 float fend2;
92 float flen1;
93 float flen2;
94 AjBool stretch;
95
96 embInit("dottup", argc, argv);
97
98 wordlen = ajAcdGetInt("wordsize");
99 seq1 = ajAcdGetSeq("bsequence");
100 seq2 = ajAcdGetSeq("asequence");
101 graph = ajAcdGetGraph("graph");
102 boxit = ajAcdGetBoolean("boxit");
103 stretch = ajAcdGetToggle("stretch");
104 xygraph = ajAcdGetGraphxy("xygraph");
105
106 begin1 = ajSeqGetBegin(seq1);
107 begin2 = ajSeqGetBegin(seq2);
108 end1 = ajSeqGetEnd(seq1);
109 end2 = ajSeqGetEnd(seq2);
110 len1 = end1 - begin1 + 1;
111 len2 = end2 - begin2 + 1;
112
113 flen1 = (float) len1;
114 flen2 = (float) len2;
115 fbegin1 = (float) begin1;
116 fbegin2 = (float) begin2;
117 fend1 = (float) end1;
118 fend2 = (float) end2;
119
120 offset1 = fbegin1;
121 offset2 = fbegin2;
122
123 ajSeqTrim(seq1);
124 ajSeqTrim(seq2);
125
126 embWordLength(wordlen);
127 if(embWordGetTable(&seq1MatchTable, seq1))
128 matchlist = embWordBuildMatchTable(seq1MatchTable, seq2, ajTrue);
129
130
131 if(stretch)
132 {
133 dottup_stretchplot(xygraph,matchlist,seq1,seq2,begin1,begin2,end1,
134 end2);
135 if(matchlist)
136 embWordMatchListDelete(&matchlist); /* free the match structures */
137 }
138
139 else
140 {
141 /* only here if stretch is false */
142
143 max= flen1;
144 if(flen2 > max)
145 max = flen2;
146
147 xmargin = ymargin = max * (float)0.15;
148
149 ajGraphOpenWin(graph, fbegin1-ymargin,fend1+ymargin,
150 fbegin2-xmargin,(float)fend2+xmargin);
151
152 ajGraphicsSetCharscale(0.5);
153
154 if(matchlist)
155 dottup_plotMatches(matchlist);
156
157 if(boxit)
158 {
159 ajGraphicsDrawposRect(fbegin1, fbegin2, fend1, fend2);
160 i = 0;
161 while(acceptableticksx[i]*numbofticks < len1)
162 i++;
163
164 if(i<=13)
165 tickgap = (float) acceptableticksx[i];
166 else
167 tickgap = (float) acceptableticksx[10];
168
169 ticklen = xmargin * (float) 0.1;
170 onefifth = xmargin * (float)0.2;
171 ajGraphicsDrawposTextAtmid(fbegin1+flen1*(float)0.5,
172 fbegin1-(onefifth*(float)3.0),
173 ajGraphGetYlabelC(graph));
174
175 if(len2/len1 > 10 )
176 {
177 /* a lot smaller then just label start and end */
178 ajGraphicsDrawposLine(fbegin1,fbegin2,fbegin1,
179 fbegin2-ticklen);
180 sprintf(ptr,"%u",ajSeqGetOffset(seq1));
181 ajGraphicsDrawposTextAtmid(fbegin1,fbegin2-(onefifth),ptr);
182
183 ajGraphicsDrawposLine(fend1,fbegin2,
184 fend1,fbegin2-ticklen);
185 sprintf(ptr,"%d",end1);
186 ajGraphicsDrawposTextAtmid(fend1,fbegin2-(onefifth),ptr);
187 }
188 else
189 for(k=fbegin1;k<fend1;k+=tickgap)
190 {
191 ajGraphicsDrawposLine(k,fbegin2,k,fbegin2-ticklen);
192 sprintf(ptr,"%d",(ajint)k);
193 ajGraphicsDrawposTextAtmid( k,fbegin2-(onefifth),ptr);
194 }
195
196 i = 0;
197 while(acceptableticks[i]*numbofticks < len2)
198 i++;
199
200 tickgap = (float) acceptableticks[i];
201 ticklen = ymargin*(float)0.1;
202 onefifth = ymargin*(float)0.2;
203 ajGraphicsDrawposTextAtlineJustify(fbegin1-(onefifth*(float)4.),
204 fbegin2+flen2*(float)0.5,
205 fbegin2-(onefifth*(float)4.),
206 fbegin2+flen2,
207 ajGraphGetXlabelC(graph),
208 0.5);
209
210 if(len1/len2 > 10 )
211 {
212 /* a lot smaller then just label start and end */
213 ajGraphicsDrawposLine(fbegin1,fbegin2,fbegin1-ticklen,
214 fbegin2);
215 sprintf(ptr,"%u",ajSeqGetOffset(seq2));
216 ajGraphicsDrawposTextAtend(fbegin1-(onefifth),fbegin2,ptr);
217
218 ajGraphicsDrawposLine(fbegin1,fend2,fbegin1-ticklen,
219 fend2);
220 sprintf(ptr,"%d",end2);
221 ajGraphicsDrawposTextAtend(fbegin2-(onefifth),fend2,ptr);
222 }
223 else
224 for(k=fbegin2;k<fend2;k+=tickgap)
225 {
226 ajGraphicsDrawposLine(fbegin1,k,fbegin1-ticklen,k);
227 sprintf(ptr,"%d",(ajint)k);
228 ajGraphicsDrawposTextAtend(fbegin1-(onefifth),k,ptr);
229 }
230 }
231 }
232
233 ajGraphicsClose();
234 ajSeqDel(&seq1);
235 ajSeqDel(&seq2);
236 ajGraphxyDel(&graph);
237 ajGraphxyDel(&xygraph);
238
239 embWordFreeTable(&seq1MatchTable);
240
241 if(matchlist)
242 embWordMatchListDelete(&matchlist); /* free the match structures */
243
244 embExit();
245
246 return 0;
247 }
248
249
250
251
252 #ifndef NO_PLOT
253 /* @funcstatic dottup_drawPlotlines *******************************************
254 **
255 ** Undocumented.
256 **
257 ** @param [r] x [void*] Undocumented
258 ** @param [r] cl [void*] Undocumented
259 ** @return [void]
260 ** @@
261 ******************************************************************************/
262
dottup_drawPlotlines(void * x,void * cl)263 static void dottup_drawPlotlines(void *x, void *cl)
264 {
265 EmbPWordMatch p;
266 PLFLT x1;
267 PLFLT y1;
268 PLFLT x2;
269 PLFLT y2;
270
271 (void) cl; /* make it used */
272
273 p = (EmbPWordMatch)x;
274
275 x1 = x2 = offset1 + (PLFLT)(p->seq1start);
276 y1 = y2 = offset2 + (PLFLT)(p->seq2start);
277
278 x2 += (PLFLT)p->length - (PLFLT)1.0;
279 y2 += (PLFLT)p->length - (PLFLT)1.0;
280
281 ajGraphicsDrawposLine(x1, y1, x2, y2);
282
283 return;
284 }
285
286
287
288
289 /* @funcstatic dottup_plotMatches *********************************************
290 **
291 ** Undocumented.
292 **
293 ** @param [r] list [const AjPList] Undocumented
294 ** @return [void]
295 ** @@
296 ******************************************************************************/
297
dottup_plotMatches(const AjPList list)298 static void dottup_plotMatches(const AjPList list)
299 {
300 ajListMapread(list, &dottup_drawPlotlines, NULL);
301
302 return;
303 }
304
305 #endif
306
307
308
309
310 /* @funcstatic dottup_stretchplot *********************************************
311 **
312 ** Undocumented.
313 **
314 ** @param [u] graph [AjPGraph] Undocumented
315 ** @param [r] matchlist [const AjPList] Undocumented
316 ** @param [r] seq1 [const AjPSeq] Undocumented
317 ** @param [r] seq2 [const AjPSeq] Undocumented
318 ** @param [r] begin1 [ajint] Undocumented
319 ** @param [r] begin2 [ajint] Undocumented
320 ** @param [r] end1 [ajint] Undocumented
321 ** @param [r] end2 [ajint] Undocumented
322 ** @return [void]
323 ** @@
324 ******************************************************************************/
325
dottup_stretchplot(AjPGraph graph,const AjPList matchlist,const AjPSeq seq1,const AjPSeq seq2,ajint begin1,ajint begin2,ajint end1,ajint end2)326 static void dottup_stretchplot(AjPGraph graph, const AjPList matchlist,
327 const AjPSeq seq1, const AjPSeq seq2,
328 ajint begin1, ajint begin2,
329 ajint end1, ajint end2)
330 {
331 EmbPWordMatch wmp = NULL;
332 float xa[1];
333 float ya[2];
334 AjPGraphdata gdata = NULL;
335 AjPStr tit = NULL;
336 float x1;
337 float y1;
338 float x2;
339 float y2;
340 AjIList iter = NULL;
341
342 tit = ajStrNew();
343 ajFmtPrintS(&tit,"%S",ajGraphGetTitleS(graph));
344
345
346 gdata = ajGraphdataNewI(1);
347 xa[0] = (float)begin1;
348 ya[0] = (float)begin2;
349
350 ajGraphSetTitleC(graph,ajStrGetPtr(tit));
351
352 ajGraphSetXlabelC(graph,ajSeqGetNameC(seq1));
353 ajGraphSetYlabelC(graph,ajSeqGetNameC(seq2));
354
355 ajGraphdataSetTypeC(gdata,"2D Plot Float");
356 ajGraphdataSetMinmax(gdata,(float)begin1,(float)end1,(float)begin2,
357 (float)end2);
358 ajGraphdataSetTruescale(gdata,(float)begin1,(float)end1,(float)begin2,
359 (float)end2);
360 ajGraphxySetXstartF(graph,(float)begin1);
361 ajGraphxySetXendF(graph,(float)end1);
362 ajGraphxySetYstartF(graph,(float)begin2);
363 ajGraphxySetYendF(graph,(float)end2);
364
365 ajGraphxySetXrangeII(graph,begin1,end1);
366 ajGraphxySetYrangeII(graph,begin2,end2);
367
368
369 if(matchlist)
370 {
371 iter = ajListIterNewread(matchlist);
372 while((wmp = ajListIterGet(iter)))
373 {
374 x1 = x2 = (float) (wmp->seq1start + begin1);
375 y1 = y2 = (float) (wmp->seq2start + begin2);
376 x2 += (float) wmp->length-1;
377 y2 += (float) wmp->length-1;
378 ajGraphAddLine(graph,x1,y1,x2,y2,0);
379 }
380 ajListIterDel(&iter);
381 }
382
383 ajGraphdataAddXY(gdata,xa,ya);
384 ajGraphDataReplace(graph,gdata);
385
386
387 ajGraphxyDisplay(graph,ajFalse);
388 ajGraphicsClose();
389
390 ajStrDel(&tit);
391
392 return;
393 }
394