Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

source: branches/items/GDE/MAFFT/mafft-7.055-with-extensions/core/io.c

Visit:

Last change on this file was 10371, checked in by aboeckma, 11 years ago
updated mafft version. Added extensions (no svn ignore, yet)
File size: 115.1 KB

Line
1	#include "mltaln.h"
2
3	static int upperCase = 0;
4
5	#define DEBUG 0
6	#define IODEBUG 0
7
8	char creverse( char f )
9	{
10	static char *table = NULL;
11	if( table == NULL )
12	{
13	int i;
14	table = AllocateCharVec(0x80);
15	for( i=0; i<0x80; i++ ) table[i] = i;
16	table['A'] = 'T';
17	table['C'] = 'G';
18	table['G'] = 'C';
19	table['T'] = 'A';
20	table['U'] = 'A';
21	table['M'] = 'K';
22	table['R'] = 'Y';
23	table['W'] = 'W';
24	table['S'] = 'S';
25	table['Y'] = 'R';
26	table['K'] = 'M';
27	table['V'] = 'B';
28	table['H'] = 'D';
29	table['D'] = 'H';
30	table['B'] = 'V';
31	table['N'] = 'N';
32	table['a'] = 't';
33	table['c'] = 'g';
34	table['g'] = 'c';
35	table['t'] = 'a';
36	table['u'] = 'a';
37	table['m'] = 'k';
38	table['r'] = 'y';
39	table['w'] = 'w';
40	table['s'] = 's';
41	table['y'] = 'r';
42	table['k'] = 'm';
43	table['v'] = 'b';
44	table['h'] = 'd';
45	table['d'] = 'h';
46	table['b'] = 'v';
47	table['n'] = 'n';
48	// table['-'] = '-';
49	// table['.'] = '.';
50	// table[''] = '';
51	}
52	return( table[(int)f] );
53	}
54
55	void sreverse( char r, char s )
56	{
57	r += strlen( s );
58	*r-- = 0;
59	while( *s )
60	r-- = creverse( s++ );
61	// r-- = ( s++ );
62	}
63
64	void gappick_samestring( char *seq )
65	{
66	char *aseq = seq;
67
68	for( ; *seq != 0; seq++ )
69	{
70	if( *seq != '-' )
71	aseq++ = seq;
72	}
73	*aseq = 0;
74	}
75
76	#if 0
77
78	static int addlocalhom2( char al1, char al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa, int skip )
79	{
80	int pos1, pos2, start1, start2, end1, end2;
81	char pt1, pt2;
82	int iscore;
83	int isumscore;
84	int sumoverlap;
85	LocalHom *tmppt;
86	int st;
87	int nlocalhom = 0;
88	pt1 = al1; pt2 = al2;
89	pos1 = off1; pos2 = off2;
90
91	isumscore = 0;
92	sumoverlap = 0;
93
94	#if 0
95	fprintf( stderr, "nlocalhom = %d in addlocalhom\n", nlocalhom );
96	fprintf( stderr, "al1 = %s, al2 = %s\n", al1, al2 );
97	fprintf( stderr, "off1 = %d, off2 = %d\n", off1, off2 );
98	fprintf( stderr, "localhopt = %p, skip = %d\n", localhompt, skip );
99	fprintf( stderr, "pt1 = \n%s\n, pt2 = \n%s\n", pt1, pt2 );
100	#endif
101
102	if( skip )
103	{
104	while( --skip > 0 ) localhompt = localhompt->next;
105	localhompt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
106	localhompt = localhompt->next;
107	// fprintf( stderr, "tmppt = %p, localhompt = %p\n", tmppt, localhompt );
108	}
109	tmppt = localhompt;
110
111	st = 0;
112	iscore = 0;
113	while( *pt1 != 0 )
114	{
115	// fprintf( stderr, "In in while loop\n" );
116	// fprintf( stderr, "pt = %c, %c, st=%d\n", pt1, pt2, st );
117	if( st == 1 && ( pt1 == '-' \|\| pt2 == '-' ) )
118	{
119	end1 = pos1 - 1;
120	end2 = pos2 - 1;
121
122	if( nlocalhom++ > 0 )
123	{
124	// fprintf( stderr, "reallocating ...\n" );
125	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
126	// fprintf( stderr, "done\n" );
127	tmppt = tmppt->next;
128	tmppt->next = NULL;
129	}
130	tmppt->start1 = start1;
131	tmppt->start2 = start2;
132	tmppt->end1 = end1 ;
133	tmppt->end2 = end2 ;
134
135	#if 1
136	isumscore += iscore;
137	sumoverlap += end2-start2+1;
138	#else
139	tmppt->overlapaa = end2-start2+1;
140	tmppt->opt = iscore * 5.8 / 600;
141	tmppt->overlapaa = overlapaa;
142	tmppt->opt = (double)opt;
143	#endif
144
145	#if 0
146	fprintf( stderr, "iscore (1)= %d\n", iscore );
147	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
148	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
149	#endif
150	iscore = 0;
151	st = 0;
152	}
153	else if( pt1 != '-' && pt2 != '-' )
154	{
155	if( st == 0 )
156	{
157	start1 = pos1; start2 = pos2;
158	st = 1;
159	}
160	iscore += n_dis[(int)amino_n[(int)pt1]][(int)amino_n[(int)pt2]];
161	// fprintf( stderr, "%c-%c, score(0) = %d\n", pt1, pt2, iscore );
162	}
163	if( *pt1++ != '-' ) pos1++;
164	if( *pt2++ != '-' ) pos2++;
165	}
166
167	if( st )
168	{
169	if( nlocalhom++ > 0 )
170	{
171	// fprintf( stderr, "reallocating ...\n" );
172	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
173	// fprintf( stderr, "done\n" );
174	tmppt = tmppt->next;
175	tmppt->next = NULL;
176	}
177	end1 = pos1 - 1;
178	end2 = pos2 - 1;
179	tmppt->start1 = start1;
180	tmppt->start2 = start2;
181	tmppt->end1 = end1 ;
182	tmppt->end2 = end2 ;
183
184	#if 1
185	isumscore += iscore;
186	sumoverlap += end2-start2+1;
187	#else
188	tmppt->overlapaa = end2-start2+1;
189	tmppt->opt = (double)iscore * 5.8 / 600;
190	tmppt->overlapaa = overlapaa;
191	tmppt->opt = (double)opt;
192	#endif
193	#if 0
194	fprintf( stderr, "score (2)= %d\n", iscore );
195	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
196	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
197	#endif
198	}
199
200	for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
201	{
202	tmppt->overlapaa = sumoverlap;
203	tmppt->opt = (double)sumscore * 5.8 / 600 / sumoverlap;
204	}
205	return( nlocalhom );
206	}
207
208	#endif
209
210
211
212	static int addlocalhom_r( char al1, char al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa, int skip )
213	{
214	int pos1, pos2, start1, start2, end1, end2;
215	char pt1, pt2;
216	double score;
217	double sumscore;
218	int sumoverlap;
219	LocalHom *tmppt = NULL; // by D.Mathog, a guess
220	int st;
221	int nlocalhom = 0;
222	pt1 = al1; pt2 = al2;
223	pos1 = off1; pos2 = off2;
224
225	sumscore = 0.0;
226	sumoverlap = 0;
227	start1 = 0; // by D.Mathog, a guess
228	start2 = 0; // by D.Mathog, a guess
229
230	#if 0
231	fprintf( stderr, "nlocalhom = %d in addlocalhom\n", nlocalhom );
232	fprintf( stderr, "al1 = %s, al2 = %s\n", al1, al2 );
233	fprintf( stderr, "off1 = %d, off2 = %d\n", off1, off2 );
234	fprintf( stderr, "localhopt = %p, skip = %d\n", localhompt, skip );
235	#endif
236	fprintf( stderr, "pt1 = \n%s\n, pt2 = \n%s\n", pt1, pt2 );
237
238	if( skip )
239	{
240	while( --skip > 0 ) localhompt = localhompt->next;
241	localhompt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
242	localhompt = localhompt->next;
243	fprintf( stderr, "tmppt = %p, localhompt = %p\n", (void )tmppt, (void )localhompt );
244	}
245	tmppt = localhompt;
246
247	st = 0;
248	score = 0.0;
249	while( *pt1 != 0 )
250	{
251	fprintf( stderr, "In in while loop\n" );
252	fprintf( stderr, "pt = %c, %c, st=%d\n", pt1, pt2, st );
253	if( st == 1 && ( pt1 == '-' \|\| pt2 == '-' ) )
254	{
255	end1 = pos1 - 1;
256	end2 = pos2 - 1;
257
258	if( nlocalhom++ > 0 )
259	{
260	// fprintf( stderr, "reallocating ...\n" );
261	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
262	// fprintf( stderr, "done\n" );
263	tmppt = tmppt->next;
264	tmppt->next = NULL;
265	}
266	tmppt->start1 = start1;
267	tmppt->start2 = start2;
268	tmppt->end1 = end1 ;
269	tmppt->end2 = end2 ;
270
271	#if 1
272	sumscore += score;
273	sumoverlap += end2-start2+1;
274	#else
275	tmppt->overlapaa = end2-start2+1;
276	tmppt->opt = score * 5.8 / 600;
277	tmppt->overlapaa = overlapaa;
278	tmppt->opt = (double)opt;
279	#endif
280
281	fprintf( stderr, "score (1)= %f\n", score );
282	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
283	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
284	score = 0.0;
285	st = 0;
286	}
287	else if( pt1 != '-' && pt2 != '-' )
288	{
289	if( st == 0 )
290	{
291	start1 = pos1; start2 = pos2;
292	st = 1;
293	}
294	score += (double)n_dis[(int)amino_n[(int)pt1]][(int)amino_n[(int)pt2]];
295	// fprintf( stderr, "%c-%c, score(0) = %f\n", pt1, pt2, score );
296	}
297	if( *pt1++ != '-' ) pos1++;
298	if( *pt2++ != '-' ) pos2++;
299	}
300	if( nlocalhom++ > 0 )
301	{
302	// fprintf( stderr, "reallocating ...\n" );
303	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
304	// fprintf( stderr, "done\n" );
305	tmppt = tmppt->next;
306	tmppt->next = NULL;
307	}
308	end1 = pos1 - 1;
309	end2 = pos2 - 1;
310	tmppt->start1 = start1;
311	tmppt->start2 = start2;
312	tmppt->end1 = end1 ;
313	tmppt->end2 = end2 ;
314
315	#if 1
316	sumscore += score;
317	sumoverlap += end2-start2+1;
318	#else
319	tmppt->overlapaa = end2-start2+1;
320	tmppt->opt = score * 5.8 / 600;
321	tmppt->overlapaa = overlapaa;
322	tmppt->opt = (double)opt;
323	#endif
324
325	fprintf( stderr, "score (2)= %f\n", score );
326	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
327	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
328
329	for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
330	{
331	tmppt->overlapaa = sumoverlap;
332	tmppt->opt = sumscore * 5.8 / 600 / sumoverlap;
333	}
334	return( nlocalhom );
335	}
336	void putlocalhom3( char al1, char al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa )
337	{
338	int pos1, pos2, start1, start2, end1, end2;
339	char pt1, pt2;
340	double score;
341	double sumscore;
342	int sumoverlap;
343	LocalHom *tmppt;
344	LocalHom *subnosento;
345	int st;
346	int saisho;
347
348	pt1 = al1; pt2 = al2;
349	pos1 = off1; pos2 = off2;
350
351	sumscore = 0.0;
352	sumoverlap = 0;
353	start1 = 0; // by Mathog, a guess
354	start2 = 0; // by Mathog, a guess
355
356	subnosento = localhompt;
357	while( subnosento->next ) subnosento = subnosento->next;
358	tmppt = subnosento;
359
360	saisho = ( localhompt->nokori == 0 );
361
362	fprintf( stderr, "localhompt = %p\n", (void *)localhompt );
363	fprintf( stderr, "tmppt = %p\n", (void *)tmppt );
364	fprintf( stderr, "subnosento = %p\n", (void *)subnosento );
365
366	st = 0;
367	score = 0.0;
368	while( *pt1 != 0 )
369	{
370	// fprintf( stderr, "pt = %c, %c, st=%d\n", pt1, pt2, st );
371	if( st == 1 && ( pt1 == '-' \|\| pt2 == '-' ) )
372	{
373	end1 = pos1 - 1;
374	end2 = pos2 - 1;
375
376	if( localhompt->nokori++ > 0 )
377	{
378	// fprintf( stderr, "reallocating ...\n" );
379	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
380	// fprintf( stderr, "done\n" );
381	tmppt = tmppt->next;
382	tmppt->next = NULL;
383	}
384	tmppt->start1 = start1;
385	tmppt->start2 = start2;
386	tmppt->end1 = end1 ;
387	tmppt->end2 = end2 ;
388
389	#if 1
390	if( divpairscore )
391	{
392	tmppt->overlapaa = end2-start2+1;
393	tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
394	}
395	else
396	{
397	sumscore += score;
398	sumoverlap += end2-start2+1;
399	}
400	#else
401	tmppt->overlapaa = overlapaa;
402	tmppt->opt = (double)opt;
403	#endif
404
405	#if 0
406	fprintf( stderr, "score (1)= %f\n", score );
407	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
408	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
409	#endif
410	score = 0.0;
411	st = 0;
412	}
413	else if( pt1 != '-' && pt2 != '-' )
414	{
415	if( st == 0 )
416	{
417	start1 = pos1; start2 = pos2;
418	st = 1;
419	}
420	score += (double)n_dis[(int)amino_n[(int)pt1]][(int)amino_n[(int)pt2]]; // - offset €Ï€€€é€Ê€€€«€â
421	// fprintf( stderr, "%c-%c, score(0) = %f\n", pt1, pt2, score );
422	}
423	if( *pt1++ != '-' ) pos1++;
424	if( *pt2++ != '-' ) pos2++;
425	}
426	if( (pt1-1) != '-' && (pt2-1) != '-' )
427	{
428	if( localhompt->nokori++ > 0 )
429	{
430	// fprintf( stderr, "reallocating ...\n" );
431	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
432	// fprintf( stderr, "done\n" );
433	tmppt = tmppt->next;
434	tmppt->next = NULL;
435	}
436
437	end1 = pos1 - 1;
438	end2 = pos2 - 1;
439	tmppt->start1 = start1;
440	tmppt->start2 = start2;
441	tmppt->end1 = end1 ;
442	tmppt->end2 = end2 ;
443
444
445	#if 1
446	if( divpairscore )
447	{
448	tmppt->overlapaa = end2-start2+1;
449	tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
450	}
451	else
452	{
453	sumscore += score;
454	sumoverlap += end2-start2+1;
455	}
456	#else
457	tmppt->overlapaa = overlapaa;
458	tmppt->opt = (double)opt;
459	#endif
460
461	#if 0
462	fprintf( stderr, "score (2)= %f\n", score );
463	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
464	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
465	#endif
466	}
467
468	fprintf( stderr, "sumscore = %f\n", sumscore );
469	if( !divpairscore )
470	{
471
472	if( !saisho ) subnosento = subnosento->next;
473	for( tmppt=subnosento; tmppt; tmppt=tmppt->next )
474	{
475	tmppt->overlapaa = sumoverlap;
476	tmppt->opt = sumscore * 5.8 / 600 / sumoverlap;
477	fprintf( stderr, "tmpptr->opt = %f\n", tmppt->opt );
478	}
479	}
480	}
481	void putlocalhom_ext( char al1, char al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa )
482	{
483	int pos1, pos2, start1, start2, end1, end2;
484	char pt1, pt2;
485	int iscore;
486	int isumscore;
487	int sumoverlap;
488	LocalHom *tmppt = localhompt;
489	int nlocalhom = 0;
490	int st;
491	pt1 = al1; pt2 = al2;
492	pos1 = off1; pos2 = off2;
493
494
495	isumscore = 0;
496	sumoverlap = 0;
497	start1 = 0; // by D.Mathog, a guess
498	start2 = 0; // by D.Mathog, a guess
499
500	st = 0;
501	iscore = 0;
502	while( *pt1 != 0 )
503	{
504	// fprintf( stderr, "pt = %c, %c, st=%d\n", pt1, pt2, st );
505	if( st == 1 && ( pt1 == '-' \|\| pt2 == '-' ) )
506	{
507	end1 = pos1 - 1;
508	end2 = pos2 - 1;
509
510	if( nlocalhom++ > 0 )
511	{
512	// fprintf( stderr, "reallocating ...\n" );
513	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
514	// fprintf( stderr, "done\n" );
515	tmppt = tmppt->next;
516	tmppt->next = NULL;
517	}
518	tmppt->start1 = start1;
519	tmppt->start2 = start2;
520	tmppt->end1 = end1 ;
521	tmppt->end2 = end2 ;
522
523	#if 1
524	if( divpairscore )
525	{
526	tmppt->overlapaa = end2-start2+1;
527	tmppt->opt = (double)iscore / tmppt->overlapaa * 5.8 / 600;
528	}
529	else
530	{
531	isumscore += iscore;
532	sumoverlap += end2-start2+1;
533	}
534	#else
535	tmppt->overlapaa = overlapaa;
536	tmppt->opt = (double)opt;
537	#endif
538
539	#if 0
540	fprintf( stderr, "iscore (1)= %d\n", iscore );
541	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
542	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
543	#endif
544	iscore = 0;
545	st = 0;
546	}
547	else if( pt1 != '-' && pt2 != '-' )
548	{
549	if( st == 0 )
550	{
551	start1 = pos1; start2 = pos2;
552	st = 1;
553	}
554	iscore += n_dis[(int)amino_n[(int)pt1]][(int)amino_n[(int)pt2]]; // - offset €Ï€€€é€Ê€€€«€â
555	// fprintf( stderr, "%c-%c, iscore(0) = %d\n", pt1, pt2, iscore );
556	}
557	if( *pt1++ != '-' ) pos1++;
558	if( *pt2++ != '-' ) pos2++;
559	}
560	if( (pt1-1) != '-' && (pt2-1) != '-' )
561	{
562	if( nlocalhom++ > 0 )
563	{
564	// fprintf( stderr, "reallocating ...\n" );
565	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
566	// fprintf( stderr, "done\n" );
567	tmppt = tmppt->next;
568	tmppt->next = NULL;
569	}
570	end1 = pos1 - 1;
571	end2 = pos2 - 1;
572	tmppt->start1 = start1;
573	tmppt->start2 = start2;
574	tmppt->end1 = end1 ;
575	tmppt->end2 = end2 ;
576
577	#if 1
578	if( divpairscore )
579	{
580	tmppt->overlapaa = end2-start2+1;
581	tmppt->opt = (double)iscore / tmppt->overlapaa * 5.8 / 600;
582	}
583	else
584	{
585	isumscore += iscore;
586	sumoverlap += end2-start2+1;
587	}
588	#else
589	tmppt->overlapaa = overlapaa;
590	tmppt->opt = (double)opt;
591	#endif
592
593	#if 0
594	fprintf( stderr, "iscore (2)= %d\n", iscore );
595	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
596	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
597	#endif
598	}
599
600	if( !divpairscore )
601	{
602	for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
603	{
604	tmppt->overlapaa = sumoverlap;
605	// tmppt->opt = (double)isumscore * 5.8 / ( 600 * sumoverlap );
606	tmppt->opt = (double)600 * 5.8 / 600;
607	// fprintf( stderr, "tmpptr->opt = %f\n", tmppt->opt );
608	}
609	}
610	}
611
612	void putlocalhom_str( char al1, char al2, double equiv, double scale, LocalHom localhompt, int off1, int off2, int opt, int overlapaa )
613	{
614	int posinaln, pos1, pos2, start1, start2, end1, end2;
615	char pt1, pt2;
616	int isumscore;
617	int sumoverlap;
618	LocalHom *tmppt = localhompt;
619	int nlocalhom = 0;
620	// int st;
621	pt1 = al1; pt2 = al2;
622	pos1 = off1; pos2 = off2;
623
624	isumscore = 0;
625	sumoverlap = 0;
626	start1 = 0; // by D.Mathog, a guess
627	start2 = 0; // by D.Mathog, a guess
628
629	posinaln = 0;
630	while( *pt1 != 0 )
631	{
632	if( pt1 != '-' && pt2 != '-' && equiv[posinaln] > 0.0 )
633	{
634	start1 = end1 = pos1; start2 = end2 = pos2;
635	if( nlocalhom++ > 0 )
636	{
637	// fprintf( stderr, "reallocating ... (posinaln=%d)\n", posinaln );
638	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
639	// fprintf( stderr, "done\n" );
640	tmppt = tmppt->next;
641	tmppt->next = NULL;
642	}
643	tmppt->start1 = start1;
644	tmppt->start2 = start2;
645	tmppt->end1 = end1 ;
646	tmppt->end2 = end2 ;
647
648	tmppt->overlapaa = 1;
649	// tmppt->opt = (double)iscore / tmppt->overlapaa * 5.8 / 600;
650	tmppt->opt = equiv[posinaln] * scale;
651	// fprintf( stdout, "pt1=%c, pt2=%c, equiv=%f\n", pt1, pt2, equiv[posinaln] );
652
653	}
654	if( *pt1++ != '-' ) pos1++;
655	if( *pt2++ != '-' ) pos2++;
656	posinaln++;
657	}
658	}
659
660
661	void putlocalhom2( char al1, char al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa )
662	{
663	int pos1, pos2, start1, start2, end1, end2;
664	char pt1, pt2;
665	int iscore;
666	int isumscore;
667	int sumoverlap;
668	LocalHom *tmppt = localhompt;
669	int nlocalhom = 0;
670	int st;
671	pt1 = al1; pt2 = al2;
672	pos1 = off1; pos2 = off2;
673
674
675	isumscore = 0;
676	sumoverlap = 0;
677	start1 = 0; // by D.Mathog, a guess
678	start2 = 0; // by D.Mathog, a guess
679
680	st = 0;
681	iscore = 0;
682	while( *pt1 != 0 )
683	{
684	// fprintf( stderr, "pt = %c, %c, st=%d\n", pt1, pt2, st );
685	if( st == 1 && ( pt1 == '-' \|\| pt2 == '-' ) )
686	{
687	end1 = pos1 - 1;
688	end2 = pos2 - 1;
689
690	if( nlocalhom++ > 0 )
691	{
692	// fprintf( stderr, "reallocating ...\n" );
693	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
694	// fprintf( stderr, "done\n" );
695	tmppt = tmppt->next;
696	tmppt->next = NULL;
697	}
698	tmppt->start1 = start1;
699	tmppt->start2 = start2;
700	tmppt->end1 = end1 ;
701	tmppt->end2 = end2 ;
702
703	#if 1
704	if( divpairscore )
705	{
706	tmppt->overlapaa = end2-start2+1;
707	tmppt->opt = (double)iscore / tmppt->overlapaa * 5.8 / 600;
708	}
709	else
710	{
711	isumscore += iscore;
712	sumoverlap += end2-start2+1;
713	}
714	#else
715	tmppt->overlapaa = overlapaa;
716	tmppt->opt = (double)opt;
717	#endif
718
719	#if 0
720	fprintf( stderr, "iscore (1)= %d\n", iscore );
721	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
722	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
723	#endif
724	iscore = 0;
725	st = 0;
726	}
727	else if( pt1 != '-' && pt2 != '-' )
728	{
729	if( st == 0 )
730	{
731	start1 = pos1; start2 = pos2;
732	st = 1;
733	}
734	iscore += n_dis[(int)amino_n[(int)pt1]][(int)amino_n[(int)pt2]]; // - offset €Ï€€€é€Ê€€€«€â
735	// fprintf( stderr, "%c-%c, iscore(0) = %d\n", pt1, pt2, iscore );
736	}
737	if( *pt1++ != '-' ) pos1++;
738	if( *pt2++ != '-' ) pos2++;
739	}
740	if( (pt1-1) != '-' && (pt2-1) != '-' )
741	{
742	if( nlocalhom++ > 0 )
743	{
744	// fprintf( stderr, "reallocating ...\n" );
745	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
746	// fprintf( stderr, "done\n" );
747	tmppt = tmppt->next;
748	tmppt->next = NULL;
749	}
750	end1 = pos1 - 1;
751	end2 = pos2 - 1;
752	tmppt->start1 = start1;
753	tmppt->start2 = start2;
754	tmppt->end1 = end1 ;
755	tmppt->end2 = end2 ;
756
757	#if 1
758	if( divpairscore )
759	{
760	tmppt->overlapaa = end2-start2+1;
761	tmppt->opt = (double)iscore / tmppt->overlapaa * 5.8 / 600;
762	}
763	else
764	{
765	isumscore += iscore;
766	sumoverlap += end2-start2+1;
767	}
768	#else
769	tmppt->overlapaa = overlapaa;
770	tmppt->opt = (double)opt;
771	#endif
772
773	#if 0
774	fprintf( stderr, "iscore (2)= %d\n", iscore );
775	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
776	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
777	#endif
778	}
779
780	if( !divpairscore )
781	{
782	for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
783	{
784	tmppt->overlapaa = sumoverlap;
785	tmppt->opt = (double)isumscore * 5.8 / ( 600 * sumoverlap );
786	// fprintf( stderr, "tmpptr->opt = %f\n", tmppt->opt );
787	}
788	}
789	}
790	void putlocalhom( char al1, char al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa )
791	{
792	int pos1, pos2, start1, start2, end1, end2;
793	char pt1, pt2;
794	double score;
795	double sumscore;
796	int sumoverlap;
797	LocalHom *tmppt = localhompt;
798	int nlocalhom = 0;
799	int st;
800	pt1 = al1; pt2 = al2;
801	pos1 = off1; pos2 = off2;
802
803
804	sumscore = 0.0;
805	sumoverlap = 0;
806	start1 = 0; // by D.Mathog, a guess
807	start2 = 0; // by D.Mathog, a guess
808
809	st = 0;
810	score = 0.0;
811	while( *pt1 != 0 )
812	{
813	// fprintf( stderr, "pt = %c, %c, st=%d\n", pt1, pt2, st );
814	if( st == 1 && ( pt1 == '-' \|\| pt2 == '-' ) )
815	{
816	end1 = pos1 - 1;
817	end2 = pos2 - 1;
818
819	if( nlocalhom++ > 0 )
820	{
821	// fprintf( stderr, "reallocating ...\n" );
822	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
823	// fprintf( stderr, "done\n" );
824	tmppt = tmppt->next;
825	tmppt->next = NULL;
826	}
827	tmppt->start1 = start1;
828	tmppt->start2 = start2;
829	tmppt->end1 = end1 ;
830	tmppt->end2 = end2 ;
831
832	#if 1
833	if( divpairscore )
834	{
835	tmppt->overlapaa = end2-start2+1;
836	tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
837	}
838	else
839	{
840	sumscore += score;
841	sumoverlap += end2-start2+1;
842	}
843	#else
844	tmppt->overlapaa = overlapaa;
845	tmppt->opt = (double)opt;
846	#endif
847
848	#if 0
849	fprintf( stderr, "score (1)= %f\n", score );
850	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
851	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
852	#endif
853	score = 0.0;
854	st = 0;
855	}
856	else if( pt1 != '-' && pt2 != '-' )
857	{
858	if( st == 0 )
859	{
860	start1 = pos1; start2 = pos2;
861	st = 1;
862	}
863	score += (double)n_dis[(int)amino_n[(int)pt1]][(int)amino_n[(int)pt2]]; // - offset €Ï€€€é€Ê€€€«€â
864	// fprintf( stderr, "%c-%c, score(0) = %f\n", pt1, pt2, score );
865	}
866	if( *pt1++ != '-' ) pos1++;
867	if( *pt2++ != '-' ) pos2++;
868	}
869	if( nlocalhom++ > 0 )
870	{
871	// fprintf( stderr, "reallocating ...\n" );
872	tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
873	// fprintf( stderr, "done\n" );
874	tmppt = tmppt->next;
875	tmppt->next = NULL;
876	}
877	end1 = pos1 - 1;
878	end2 = pos2 - 1;
879	tmppt->start1 = start1;
880	tmppt->start2 = start2;
881	tmppt->end1 = end1 ;
882	tmppt->end2 = end2 ;
883
884	#if 1
885	if( divpairscore )
886	{
887	tmppt->overlapaa = end2-start2+1;
888	tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
889	}
890	else
891	{
892	sumscore += score;
893	sumoverlap += end2-start2+1;
894	}
895	#else
896	tmppt->overlapaa = overlapaa;
897	tmppt->opt = (double)opt;
898	#endif
899
900	#if 0
901	fprintf( stderr, "score (2)= %f\n", score );
902	fprintf( stderr, "al1: %d - %d\n", start1, end1 );
903	fprintf( stderr, "al2: %d - %d\n", start2, end2 );
904	#endif
905
906	if( !divpairscore )
907	{
908	for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
909	{
910	tmppt->overlapaa = sumoverlap;
911	tmppt->opt = sumscore * 5.8 / 600 / sumoverlap;
912	// fprintf( stderr, "tmpptr->opt = %f\n", tmppt->opt );
913	}
914	}
915	}
916
917	char cutal( char al, int al_display_start, int start, int end )
918	{
919	int pos;
920	char *pt = al;
921	char *val = NULL;
922
923	pos = al_display_start;
924	do
925	{
926	if( start == pos ) val = pt;
927	if( end == pos ) break;
928	// fprintf( stderr, "pos=%d, pt=%c, val=%p\n", pos, pt, val );
929	if( *pt != '-' ) pos++;
930	} while( *pt++ != 0 );
931	*(pt+1) = 0;
932	return( val );
933	}
934
935	void ErrorExit( char *message )
936	{
937	fprintf( stderr, "%s\n", message );
938	exit( 1 );
939	}
940
941	void strncpy_caseC( char str1, char str2, int len )
942	{
943	if( dorp == 'd' && upperCase > 0 )
944	{
945	while( len-- )
946	str1++ = toupper( str2++ );
947	}
948	else strncpy( str1, str2, len );
949	}
950
951	void seqUpper( int nseq, char *seq ) / not used */
952	{
953	int i, j, len;
954	for( i=0; i<nseq; i++ )
955	{
956	len = strlen( seq[i] );
957	for( j=0; j<len; j++ )
958	seq[i][j] = toupper( seq[i][j] );
959	}
960	}
961
962	void seqLower( int nseq, char **seq )
963	{
964	int i, j, len;
965	for( i=0; i<nseq; i++ )
966	{
967	len = strlen( seq[i] );
968	for( j=0; j<len; j++ )
969	seq[i][j] = tolower( seq[i][j] );
970	}
971	}
972
973	int getaline_fp_eof( char s, int l, FILE fp ) /* end of file -> return 1 */
974	{
975	int c, i = 0 ;
976	int noteofflag = 0;
977	for( i=0; i<l && ( noteofflag = ( (c=getc(fp)) != EOF ) ) && c != '\n'; i++ )
978	*s++ = c;
979	*s = '\0' ;
980	return( !noteofflag );
981	}
982
983	int getaline_fp_eof_new(s, l, fp) /* end of file -> return 1 */
984	char s[] ; int l ; FILE *fp ;
985	{
986	int c = 0, i = 0 ;
987	int noteofflag = 0;
988
989	if( feof( fp ) ) return( 1 );
990
991	for( i=0; i<l && ( noteofflag = ( (c=getc(fp)) != EOF ) ) && c != '\n'; i++ )
992	{ *s++ = c ; }
993	*s = '\0' ;
994	if( c != '\n' && c != EOF ) while( getc(fp) != '\n' )
995	;
996	return( !noteofflag );
997	}
998
999	int myfgets(s, l, fp) /* l°ÊŸå€Ï¡¢¹ÔËö€Þ€ÇÆÉ€ßÈô€Ð€¹ */
1000	char s[] ; int l ; FILE *fp ;
1001	{
1002	int c = 0, i = 0 ;
1003
1004	if( feof( fp ) ) return( 1 );
1005
1006	for( i=0; i<l && ( c=getc( fp ) ) != '\n'; i++ )
1007	*s++ = c;
1008	*s = '\0' ;
1009	if( c != '\n' )
1010	while( getc(fp) != '\n' )
1011	;
1012	return( 0 );
1013	}
1014
1015	float input_new( FILE *fp, int d )
1016	{
1017	char mojiretsu[10];
1018	int i, c;
1019
1020	c = getc( fp );
1021	if( c != '\n' )
1022	ungetc( c, fp );
1023
1024	for( i=0; i<d; i++ )
1025	mojiretsu[i] = getc( fp );
1026	mojiretsu[i] = 0;
1027
1028	return( atof( mojiretsu ) );
1029	}
1030
1031
1032	void PreRead( FILE fp, int locnjob, int *locnlenmax )
1033	{
1034	int i, nleni;
1035	char b[B];
1036
1037	fgets( b, B-1, fp ); *locnjob = atoi( b );
1038	*locnlenmax = 0;
1039	i=0;
1040	while( i<*locnjob )
1041	{
1042	fgets( b, B-1, fp );
1043	if( !strncmp( b, "=", 1 ) )
1044	{
1045	i++;
1046	fgets( b, B-1, fp ); nleni = atoi( b );
1047	if( nleni > locnlenmax ) locnlenmax = nleni;
1048	}
1049	}
1050	if( *locnlenmax > N )
1051	{
1052	fprintf( stderr, "TOO LONG SEQUENCE!\n" );
1053	exit( 1 );
1054	}
1055	if( njob > M )
1056	{
1057	fprintf( stderr, "TOO MANY SEQUENCE!\n" );
1058	fprintf( stderr, "%d > %d\n", njob, M );
1059	exit( 1 );
1060	}
1061	}
1062
1063	int allSpace( char *str )
1064	{
1065	int value = 1;
1066	while( str ) value = ( !isdigit( *str++ ) );
1067	return( value );
1068	}
1069
1070	void Read( char name[M][B], int nlen[M], char **seq )
1071	{
1072	extern void FRead( FILE x, char y[M][B], int z[M], char *w );
1073	FRead( stdin, name, nlen, seq );
1074	}
1075
1076
1077	void FRead( FILE fp, char name[][B], int nlen[], char *seq )
1078	{
1079	int i, j;
1080	char b[B];
1081
1082	fgets( b, B-1, fp );
1083	#if DEBUG
1084	fprintf( stderr, "b = %s\n", b );
1085	#endif
1086
1087	if( strstr( b, "onnet" ) ) scoremtx = 1;
1088	else if( strstr( b, "DnA" ) )
1089	{
1090	scoremtx = -1;
1091	upperCase = -1;
1092	}
1093	else if( strstr( b, "dna" ) )
1094	{
1095	scoremtx = -1;
1096	upperCase = 0;
1097	}
1098	else if( strstr( b, "DNA" ) )
1099	{
1100	scoremtx = -1;
1101	upperCase = 1;
1102	}
1103	else if( strstr( b, "M-Y" ) \|\| strstr( b, "iyata" ) ) scoremtx = 2;
1104	else scoremtx = 0;
1105	#if DEBUG
1106	fprintf( stderr, " %s->scoremtx = %d\n", b, scoremtx );
1107	#endif
1108
1109	geta2 = GETA2;
1110
1111	#if 0
1112	if( strlen( b ) >=25 )
1113	{
1114	b[25] = 0;
1115	#if DEBUG
1116	fprintf( stderr, "kimuraR = %s\n", b+20 );
1117	#endif
1118	kimuraR = atoi( b+20 );
1119
1120	if( kimuraR < 0 \|\| 20 < kimuraR ) ErrorExit( "Illeagal kimuraR value.\n" );
1121	if( allSpace( b+20 ) ) kimuraR = NOTSPECIFIED;
1122	}
1123	else kimuraR = NOTSPECIFIED;
1124	#if DEBUG
1125	fprintf( stderr, "kimuraR = %d\n", kimuraR );
1126	#endif
1127
1128	if( strlen( b ) >=20 )
1129	{
1130	b[20] = 0;
1131	#if DEBUG
1132	fprintf( stderr, "pamN = %s\n", b+15 );
1133	#endif
1134	pamN = atoi( b+15 );
1135	if( pamN < 0 \|\| 400 < pamN ) ErrorExit( "Illeagal pam value.\n" );
1136	if( allSpace( b+15 ) ) pamN = NOTSPECIFIED;
1137	}
1138	else pamN = NOTSPECIFIED;
1139
1140	if( strlen( b ) >= 15 )
1141	{
1142	b[15] = 0;
1143	#if DEBUG
1144	fprintf( stderr, "poffset = %s\n", b+10 );
1145	#endif
1146	poffset = atoi( b+10 );
1147	if( poffset > 500 ) ErrorExit( "Illegal extending gap ppenalty\n" );
1148	if( allSpace( b+10 ) ) poffset = NOTSPECIFIED;
1149	}
1150	else poffset = NOTSPECIFIED;
1151
1152	if( strlen( b ) >= 10 )
1153	{
1154	b[10] = 0;
1155	#if DEBUG
1156	fprintf( stderr, "ppenalty = %s\n", b+5 );
1157	#endif
1158	ppenalty = atoi( b+5 );
1159	if( ppenalty > 0 ) ErrorExit( "Illegal opening gap ppenalty\n" );
1160	if( allSpace( b+5 ) ) ppenalty = NOTSPECIFIED;
1161	}
1162	else ppenalty = NOTSPECIFIED;
1163	#endif
1164
1165	for( i=0; i<njob; i++ )
1166	{
1167	getaline_fp_eof_new( b, B-1, fp );
1168	strcpy( name[i], b );
1169	#if DEBUG
1170	fprintf( stderr, "name[%d] = %s\n", i, name[i] );
1171	#endif
1172	fgets( b, B-1, fp ); nlen[i] = atoi( b ); /* seq i no nagasa */
1173	seq[i][0] = 0;
1174	if( nlen[i] ) for( j=0; j <= (nlen[i]-1)/C; j++ )
1175	{
1176	getaline_fp_eof_new( b, B-1, fp );
1177	/* b[C] = 0; */
1178	strcat( seq[i], b );
1179	}
1180	seq[i][nlen[i]] = 0;
1181	}
1182	if( scoremtx == -1 && upperCase != -1 ) seqLower( njob, seq );
1183	}
1184
1185
1186	static int countKUorWA( FILE *fp )
1187	{
1188	int value;
1189	int c, b;
1190
1191	value= 0;
1192	b = '\n';
1193	while( ( c = getc( fp ) ) != EOF )
1194	{
1195	if( b == '\n' && ( c == '>' ) )
1196	value++;
1197	b = c;
1198	}
1199	rewind( fp );
1200	return( value );
1201	}
1202
1203	void searchKUorWA( FILE *fp )
1204	{
1205	int c, b;
1206	b = '\n';
1207	while( !( ( ( c = getc( fp ) ) == '>' \|\| c == EOF ) && b == '\n' ) )
1208	b = c;
1209	ungetc( c, fp );
1210	}
1211
1212	static int onlyGraph( char *str )
1213	{
1214	char tmp;
1215	char *res = str;
1216	char *bk = str;
1217
1218	// while( (tmp=str++) ) if( isgraph( tmp ) ) res++ = tmp;
1219	while( (tmp=*str++) )
1220	{
1221	if( 0x20 < tmp && tmp < 0x7f ) *res++ = tmp;
1222	if( tmp == '>' )
1223	{
1224	fprintf( stderr, "========================================================\n" );
1225	fprintf( stderr, "========================================================\n" );
1226	fprintf( stderr, "=== \n" );
1227	fprintf( stderr, "=== ERROR!! \n" );
1228	fprintf( stderr, "=== In the '--anysymbol' and '--preservecase' modes, \n" );
1229	fprintf( stderr, "=== '>' in sequence is unacceptable.\n" );
1230	fprintf( stderr, "=== \n" );
1231	fprintf( stderr, "========================================================\n" );
1232	fprintf( stderr, "========================================================\n" );
1233	exit( 1 );
1234	}
1235	}
1236	*res = 0;
1237	return( res - bk );
1238	}
1239
1240	static int onlyAlpha_lower( char *str )
1241	{
1242	char tmp;
1243	char *res = str;
1244	char *bk = str;
1245
1246	while( (tmp=*str++) )
1247	if( isalpha( tmp ) \|\| tmp == '-' \|\| tmp == '*' \|\| tmp == '.' )
1248	*res++ = tolower( tmp );
1249	*res = 0;
1250	return( res - bk );
1251	}
1252	static int onlyAlpha_upper( char *str )
1253	{
1254	char tmp;
1255	char *res = str;
1256	char *bk = str;
1257
1258	while( (tmp=*str++) )
1259	if( isalpha( tmp ) \|\| tmp == '-' \|\| tmp == '*' \|\| tmp == '.' )
1260	*res++ = toupper( tmp );
1261	*res = 0;
1262	return( res - bk );
1263	}
1264
1265	void kake2hiku( char *str )
1266	{
1267	do
1268	if( str == '' ) *str = '-';
1269	while( *str++ );
1270	}
1271
1272	char load1SeqWithoutName_realloc_casepreserve( FILE fpp )
1273	{
1274	int c, b;
1275	char *cbuf;
1276	int size = N;
1277	char *val;
1278
1279	val = malloc( (size+1) * sizeof( char ) );
1280	cbuf = val;
1281
1282	b = '\n';
1283	while( ( c = getc( fpp ) ) != EOF &&
1284	!( ( c == '>' \|\| c == '(' \|\| c == EOF ) && b == '\n' ) )
1285	{
1286	cbuf++ = (char)c; / Ä¹€¹€®€Æ€â€·€é€Ê€€ */
1287	if( cbuf - val == size )
1288	{
1289	size += N;
1290	fprintf( stderr, "reallocating...\n" );
1291	val = (char )realloc( val, (size+1) sizeof( char ) );
1292	if( !val )
1293	{
1294	fprintf( stderr, "Allocation error in load1SeqWithoutName_realloc \n" );
1295	exit( 1 );
1296	}
1297	fprintf( stderr, "done.\n" );
1298	cbuf = val + size-N;
1299	}
1300	b = c;
1301	}
1302	ungetc( c, fpp );
1303	*cbuf = 0;
1304	onlyGraph( val );
1305	// kake2hiku( val );
1306	return( val );
1307	}
1308
1309	char load1SeqWithoutName_realloc( FILE fpp )
1310	{
1311	int c, b;
1312	char *cbuf;
1313	int size = N;
1314	char *val;
1315
1316	val = malloc( (size+1) * sizeof( char ) );
1317	cbuf = val;
1318
1319	b = '\n';
1320	while( ( c = getc( fpp ) ) != EOF &&
1321	!( ( c == '>' \|\| c == '(' \|\| c == EOF ) && b == '\n' ) )
1322	{
1323	cbuf++ = (char)c; / Ä¹€¹€®€Æ€â€·€é€Ê€€ */
1324	if( cbuf - val == size )
1325	{
1326	size += N;
1327	fprintf( stderr, "reallocating...\n" );
1328	val = (char )realloc( val, (size+1) sizeof( char ) );
1329	if( !val )
1330	{
1331	fprintf( stderr, "Allocation error in load1SeqWithoutName_realloc \n" );
1332	exit( 1 );
1333	}
1334	fprintf( stderr, "done.\n" );
1335	cbuf = val + size-N;
1336	}
1337	b = c;
1338	}
1339	ungetc( c, fpp );
1340	*cbuf = 0;
1341	if( dorp == 'd' )
1342	onlyAlpha_lower( val );
1343	else
1344	onlyAlpha_upper( val );
1345	kake2hiku( val );
1346	return( val );
1347	}
1348
1349	int load1SeqWithoutName_new( FILE fpp, char cbuf )
1350	{
1351	int c, b;
1352	char *bk = cbuf;
1353
1354	b = '\n';
1355	while( ( c = getc( fpp ) ) != EOF && /* by T. Nishiyama */
1356	!( ( c == '>' \|\| c == '(' \|\| c == EOF ) && b == '\n' ) )
1357	{
1358	cbuf++ = (char)c; / Ä¹€¹€®€Æ€â€·€é€Ê€€ */
1359	b = c;
1360	}
1361	ungetc( c, fpp );
1362	*cbuf = 0;
1363	if( dorp == 'd' )
1364	onlyAlpha_lower( bk );
1365	else
1366	onlyAlpha_upper( bk );
1367	kake2hiku( bk );
1368	return( 0 );
1369	}
1370
1371
1372	void readDataforgaln( FILE fp, char name, int nlen, char **seq )
1373	{
1374	int i;
1375	static char *tmpseq = NULL;
1376
1377	#if 0
1378	if( !tmpseq )
1379	{
1380	tmpseq = AllocateCharVec( N );
1381	}
1382	#endif
1383
1384	rewind( fp );
1385	searchKUorWA( fp );
1386
1387	for( i=0; i<njob; i++ )
1388	{
1389	name[i][0] = '='; getc( fp );
1390	#if 0
1391	fgets( name[i]+1, B-2, fp );
1392	j = strlen( name[i] );
1393	if( name[i][j-1] != '\n' )
1394	ErrorExit( "Too long name\n" );
1395	name[i][j-1] = 0;
1396	#else
1397	myfgets( name[i]+1, B-2, fp );
1398	#endif
1399	#if 0
1400	fprintf( stderr, "name[%d] = %s\n", i, name[i] );
1401	#endif
1402	tmpseq = load1SeqWithoutName_realloc( fp );
1403	strcpy( seq[i], tmpseq );
1404	nlen[i] = strlen( seq[i] );
1405	free( tmpseq );
1406	}
1407	if( dorp == 'd' && upperCase != -1 ) seqLower( njob, seq );
1408	#if 0
1409	free( tmpseq );
1410	#endif
1411	}
1412
1413	void readData_varlen( FILE fp, char name, int nlen, char **seq )
1414	{
1415	int i;
1416	static char *tmpseq = NULL;
1417
1418	rewind( fp );
1419	searchKUorWA( fp );
1420
1421	for( i=0; i<njob; i++ )
1422	{
1423	name[i][0] = '='; getc( fp );
1424	#if 0
1425	fgets( name[i]+1, B-2, fp );
1426	j = strlen( name[i] );
1427	if( name[i][j-1] != '\n' )
1428	ErrorExit( "Too long name\n" );
1429	name[i][j-1] = 0;
1430	#else
1431	myfgets( name[i]+1, B-2, fp );
1432	#endif
1433	#if 0
1434	fprintf( stderr, "name[%d] = %s\n", i, name[i] );
1435	#endif
1436	tmpseq = load1SeqWithoutName_realloc( fp );
1437	nlen[i] = strlen( tmpseq );
1438	// fprintf( stderr, "nlen[%d] = %d\n", i+1, nlen[i] );
1439	seq[i] = calloc( nlen[i]+1, sizeof( char ) );
1440	strcpy( seq[i], tmpseq );
1441	free( tmpseq );
1442	}
1443	if( dorp == 'd' && upperCase != -1 ) seqLower( njob, seq );
1444	#if 0
1445	free( tmpseq );
1446	#endif
1447	}
1448
1449	void readData_pointer2( FILE fp, int nseq, char name, int nlen, char **seq )
1450	{
1451	int i;
1452	static char *tmpseq = NULL;
1453
1454	#if 0
1455	if( !tmpseq )
1456	{
1457	tmpseq = AllocateCharVec( N );
1458	}
1459	#endif
1460
1461	rewind( fp );
1462	searchKUorWA( fp );
1463
1464	for( i=0; i<nseq; i++ )
1465	{
1466	name[i][0] = '='; getc( fp );
1467	#if 0
1468	fgets( name[i]+1, B-2, fp );
1469	j = strlen( name[i] );
1470	if( name[i][j-1] != '\n' )
1471	ErrorExit( "Too long name\n" );
1472	name[i][j-1] = 0;
1473	#else
1474	myfgets( name[i]+1, B-2, fp );
1475	#endif
1476	#if 0
1477	fprintf( stderr, "name[%d] = %s\n", i, name[i] );
1478	#endif
1479	tmpseq = load1SeqWithoutName_realloc( fp );
1480	strcpy( seq[i], tmpseq );
1481	free( tmpseq );
1482	nlen[i] = strlen( seq[i] );
1483	}
1484	if( dorp == 'd' && upperCase != -1 ) seqLower( nseq, seq );
1485	#if 0
1486	free( tmpseq );
1487	#endif
1488	if( outnumber )
1489	{
1490	char *namebuf;
1491	char *cptr;
1492	namebuf = calloc( B+100, sizeof( char ) );
1493	for( i=0; i<nseq; i++ )
1494	{
1495	namebuf[0] = '=';
1496	cptr = strstr( name[i], "_numo_e_" );
1497	if( cptr )
1498	sprintf( namebuf+1, "_numo_s_%08d_numo_e_%s", i+1, cptr+8 );
1499	else
1500	sprintf( namebuf+1, "_numo_s_%08d_numo_e_%s", i+1, name[i]+1 );
1501	strncpy( name[i], namebuf, B );
1502	name[i][B-1] = 0;
1503	}
1504	free( namebuf );
1505	// exit( 1 );
1506	}
1507	}
1508
1509
1510	void readData_pointer_casepreserve( FILE fp, char name, int nlen, char **seq )
1511	{
1512	int i;
1513	static char *tmpseq = NULL;
1514
1515	#if 0
1516	if( !tmpseq )
1517	{
1518	tmpseq = AllocateCharVec( N );
1519	}
1520	#endif
1521
1522	rewind( fp );
1523	searchKUorWA( fp );
1524
1525	for( i=0; i<njob; i++ )
1526	{
1527	name[i][0] = '='; getc( fp );
1528	#if 0
1529	fgets( name[i]+1, B-2, fp );
1530	j = strlen( name[i] );
1531	if( name[i][j-1] != '\n' )
1532	ErrorExit( "Too long name\n" );
1533	name[i][j-1] = 0;
1534	#else
1535	myfgets( name[i]+1, B-2, fp );
1536	#endif
1537	#if 0
1538	fprintf( stderr, "name[%d] = %s\n", i, name[i] );
1539	#endif
1540	tmpseq = load1SeqWithoutName_realloc_casepreserve( fp );
1541	strcpy( seq[i], tmpseq );
1542	free( tmpseq );
1543	nlen[i] = strlen( seq[i] );
1544	}
1545	}
1546
1547
1548	void readData_pointer( FILE fp, char name, int nlen, char **seq )
1549	{
1550	int i;
1551	static char *tmpseq = NULL;
1552
1553	#if 0
1554	if( !tmpseq )
1555	{
1556	tmpseq = AllocateCharVec( N );
1557	}
1558	#endif
1559
1560	rewind( fp );
1561	searchKUorWA( fp );
1562
1563	for( i=0; i<njob; i++ )
1564	{
1565	name[i][0] = '='; getc( fp );
1566	#if 0
1567	fgets( name[i]+1, B-2, fp );
1568	j = strlen( name[i] );
1569	if( name[i][j-1] != '\n' )
1570	ErrorExit( "Too long name\n" );
1571	name[i][j-1] = 0;
1572	#else
1573	myfgets( name[i]+1, B-2, fp );
1574	#endif
1575	#if 0
1576	fprintf( stderr, "name[%d] = %s\n", i, name[i] );
1577	#endif
1578	tmpseq = load1SeqWithoutName_realloc( fp );
1579	strcpy( seq[i], tmpseq );
1580	free( tmpseq );
1581	nlen[i] = strlen( seq[i] );
1582	}
1583	if( dorp == 'd' && upperCase != -1 ) seqLower( njob, seq );
1584	#if 0
1585	free( tmpseq );
1586	#endif
1587	if( outnumber )
1588	{
1589	char *namebuf;
1590	char *cptr;
1591	namebuf = calloc( B+100, sizeof( char ) );
1592	for( i=0; i<njob; i++ )
1593	{
1594	namebuf[0] = '=';
1595	cptr = strstr( name[i], "_numo_e_" );
1596	if( cptr )
1597	sprintf( namebuf+1, "_numo_s_%08d_numo_e_%s", i+1, cptr+8 );
1598	else
1599	sprintf( namebuf+1, "_numo_s_%08d_numo_e_%s", i+1, name[i]+1 );
1600	strncpy( name[i], namebuf, B );
1601	name[i][B-1] = 0;
1602	}
1603	free( namebuf );
1604	// exit( 1 );
1605	}
1606	}
1607
1608	void readData( FILE fp, char name[][B], int nlen[], char *seq )
1609	{
1610	int i;
1611	static char *tmpseq = NULL;
1612
1613	#if 0
1614	if( !tmpseq )
1615	{
1616	tmpseq = AllocateCharVec( N );
1617	}
1618	#endif
1619
1620	rewind( fp );
1621	searchKUorWA( fp );
1622
1623	for( i=0; i<njob; i++ )
1624	{
1625	name[i][0] = '='; getc( fp );
1626	#if 0
1627	fgets( name[i]+1, B-2, fp );
1628	j = strlen( name[i] );
1629	if( name[i][j-1] != '\n' )
1630	ErrorExit( "Too long name\n" );
1631	name[i][j-1] = 0;
1632	#else
1633	myfgets( name[i]+1, B-2, fp );
1634	#endif
1635	#if 0
1636	fprintf( stderr, "name[%d] = %s\n", i, name[i] );
1637	#endif
1638	tmpseq = load1SeqWithoutName_realloc( fp );
1639	strcpy( seq[i], tmpseq );
1640	nlen[i] = strlen( seq[i] );
1641	free( tmpseq );
1642	}
1643	if( dorp == 'd' && upperCase != -1 ) seqLower( njob, seq );
1644	#if 0
1645	free( tmpseq );
1646	#endif
1647	}
1648
1649	void cutAlignment( FILE fp, int regtable, char revtable, int outtable, char name, char outseq )
1650	{
1651	int i, j;
1652	int outlen;
1653	static char *tmpseq = NULL;
1654	static char *dumname = NULL;
1655	char fs, rs;
1656	int npos, lpos;
1657	int startpos, endpos, seqlen;
1658
1659	if( dumname == NULL )
1660	{
1661	dumname = AllocateCharVec( N );
1662	}
1663
1664	rewind( fp );
1665	searchKUorWA( fp );
1666
1667
1668	npos = 0;
1669	for( i=0; i<njob; i++ )
1670	{
1671	dumname[0] = '>'; getc( fp );
1672	myfgets( dumname+1, B-1, fp );
1673	tmpseq = load1SeqWithoutName_realloc_casepreserve( fp );
1674
1675	if( outtable[i] )
1676	{
1677	// putc( '>', stdout );
1678	// puts( dumname+1 );
1679
1680
1681	strncat( name[npos], dumname, B-1 );
1682	name[npos][B-1] = 0;
1683
1684	if( dorp == 'd' && upperCase != -1 ) seqLower( 1, &tmpseq );
1685	seqlen = strlen( tmpseq );
1686	lpos = 0;
1687	for( j=0; j<5; j++ )
1688	{
1689	if( regtable[0][j2] == -1 && regtable[0][j2+1] == -1 ) continue;
1690
1691	startpos = regtable[0][j*2];
1692	endpos = regtable[0][j*2+1];
1693	if( startpos > endpos )
1694	{
1695	endpos = regtable[0][j*2];
1696	startpos = regtable[0][j*2+1];
1697	}
1698
1699	if( startpos < 0 ) startpos = 0;
1700	if( endpos < 0 ) endpos = 0;
1701	if( endpos >= seqlen ) endpos = seqlen-1;
1702	if( startpos >= seqlen ) startpos = seqlen-1;
1703
1704	// fprintf( stderr, "startpos = %d, endpos = %d\n", startpos, endpos );
1705
1706	outlen = endpos - startpos+1;
1707	if( revtable[0][j] == 'f' )
1708	{
1709	// fprintf( stderr, "regtable[%d][st] = %d\n", i, regtable[0][j*2+0] );
1710	// fprintf( stderr, "regtable[%d][en] = %d\n", i, regtable[0][j*2+1] );
1711	// fprintf( stderr, "outlen = %d\n", outlen );
1712	// fprintf( stdout, "%.s\n", outlen, tmpseq+regtable[0][j2] );
1713	strncpy( outseq[npos] + lpos, tmpseq+startpos, outlen );
1714	lpos += outlen;
1715	}
1716	else
1717	{
1718	fs = AllocateCharVec( outlen+1 );
1719	rs = AllocateCharVec( outlen+1 );
1720
1721	fs[outlen] = 0;
1722	strncpy( fs, tmpseq+startpos, outlen );
1723	sreverse( rs, fs );
1724	// fprintf( stdout, "%s\n", rs );
1725	strncpy( outseq[npos] + lpos, rs, outlen );
1726	lpos += outlen;
1727	free( fs );
1728	free( rs );
1729	}
1730	outseq[npos][lpos] = 0;
1731	}
1732	npos++;
1733	}
1734	free( tmpseq );
1735	}
1736	}
1737
1738	void cutData( FILE fp, int regtable, char revtable, int outtable )
1739	{
1740	int i, j;
1741	int outlen, seqlen, startpos, endpos;
1742	static char *tmpseq = NULL;
1743	static char *dumname = NULL;
1744	char fs, rs;
1745
1746	if( dumname == NULL )
1747	{
1748	dumname = AllocateCharVec( N );
1749	}
1750
1751	rewind( fp );
1752	searchKUorWA( fp );
1753
1754	for( i=0; i<njob; i++ )
1755	{
1756	dumname[0] = '='; getc( fp );
1757	myfgets( dumname+1, B-2, fp );
1758	tmpseq = load1SeqWithoutName_realloc_casepreserve( fp );
1759
1760	if( outtable[i] )
1761	{
1762	gappick_samestring( tmpseq );
1763	putc( '>', stdout );
1764	puts( dumname+1 );
1765
1766	seqlen = strlen( tmpseq );
1767
1768	if( dorp == 'd' && upperCase != -1 ) seqLower( 1, &tmpseq );
1769	if( outtable[i] == 2 )
1770	{
1771	startpos = 0;
1772	endpos = seqlen-1;
1773	outlen = endpos - startpos + 1;
1774	fprintf( stdout, "%.*s\n", outlen, tmpseq+startpos );
1775	}
1776	else
1777	{
1778	for( j=0; j<5; j++ )
1779	{
1780	if( regtable[i][j2] == -1 && regtable[i][j2+1] == -1 ) continue;
1781
1782	startpos = regtable[i][j*2];
1783	endpos = regtable[i][j*2+1];
1784
1785	if( startpos > endpos )
1786	{
1787	endpos = regtable[i][j*2];
1788	startpos = regtable[i][j*2+1];
1789	}
1790
1791	if( startpos < 0 ) startpos = 0;
1792	if( endpos < 0 ) endpos = 0;
1793	if( endpos >= seqlen ) endpos = seqlen-1;
1794	if( startpos >= seqlen ) startpos = seqlen-1;
1795
1796	outlen = endpos - startpos + 1;
1797	if( revtable[i][j] == 'f' )
1798	{
1799	fprintf( stderr, "startpos = %d\n", startpos );
1800	fprintf( stderr, "endpos = %d\n", endpos );
1801	fprintf( stderr, "outlen = %d\n", outlen );
1802	fprintf( stdout, "%.*s\n", outlen, tmpseq+startpos );
1803	}
1804	else
1805	{
1806	fs = AllocateCharVec( outlen+1 );
1807	rs = AllocateCharVec( outlen+1 );
1808
1809	fs[outlen] = 0;
1810	strncpy( fs, tmpseq+startpos, outlen );
1811	sreverse( rs, fs );
1812	fprintf( stdout, "%s\n", rs );
1813	free( fs );
1814	free( rs );
1815	}
1816	}
1817	}
1818	}
1819	free( tmpseq );
1820	}
1821	}
1822
1823	void catData( FILE *fp )
1824	{
1825	int i;
1826	static char *tmpseq = NULL;
1827	static char *dumname = NULL;
1828	// char *cptr;
1829
1830	if( dumname == NULL )
1831	{
1832	dumname = AllocateCharVec( N );
1833	}
1834
1835	rewind( fp );
1836	searchKUorWA( fp );
1837
1838	for( i=0; i<njob; i++ )
1839	{
1840	dumname[0] = '='; getc( fp );
1841	myfgets( dumname+1, B-2, fp );
1842	if( outnumber )
1843	{
1844	fprintf( stdout, ">_numo_s_%08d_numo_e_", i+1 );
1845	}
1846	else
1847	{
1848	putc( '>', stdout );
1849	}
1850	puts( dumname+1 );
1851	tmpseq = load1SeqWithoutName_realloc( fp );
1852	if( dorp == 'd' && upperCase != -1 ) seqLower( 1, &tmpseq );
1853	puts( tmpseq );
1854	free( tmpseq );
1855	}
1856	}
1857
1858	int countATGC( char s, int total )
1859	{
1860	int nATGC;
1861	int nChar;
1862	char c;
1863	nATGC = nChar = 0;
1864
1865	if( *s == 0 )
1866	{
1867	*total = 0;
1868	return( 0 );
1869	}
1870
1871	do
1872	{
1873	c = tolower( *s );
1874	if( isalpha( c ) )
1875	{
1876	nChar++;
1877	if( c == 'a' \|\| c == 't' \|\| c == 'g' \|\| c == 'c' \|\| c == 'u' \|\| c == 'n' )
1878	nATGC++;
1879	}
1880	}
1881	while( *++s );
1882
1883	*total = nChar;
1884	return( nATGC );
1885	}
1886
1887	double countATGCbk( char *s )
1888	{
1889	int nATGC;
1890	int nChar;
1891	char c;
1892	nATGC = nChar = 0;
1893
1894	do
1895	{
1896	c = tolower( *s );
1897	if( isalpha( c ) )
1898	{
1899	nChar++;
1900	if( c == 'a' \|\| c == 't' \|\| c == 'g' \|\| c == 'c' \|\| c == 'u' \|\| c == 'n' )
1901	nATGC++;
1902	}
1903	}
1904	while( *++s );
1905	return( (double)nATGC / nChar );
1906	}
1907
1908
1909	int countnogaplen( char *seq )
1910	{
1911	int val = 0;
1912	while( *seq )
1913	if( *seq++ != '-' ) val++;
1914	return( val );
1915	}
1916
1917	int countnormalletters( char seq, char ref )
1918	{
1919	int val = 0;
1920	while( *seq )
1921	if( strchr( ref, *seq++ ) ) val++;
1922	return( val );
1923	}
1924
1925	void getnumlen_casepreserve( FILE fp, int nlenminpt )
1926	{
1927	int total;
1928	int nsite = 0;
1929	int atgcnum;
1930	int i, tmp;
1931	char tmpseq, tmpname;
1932	double atgcfreq;
1933	tmpname = AllocateCharVec( N );
1934	njob = countKUorWA( fp );
1935	searchKUorWA( fp );
1936	nlenmax = 0;
1937	*nlenminpt = 99999999;
1938	atgcnum = 0;
1939	total = 0;
1940	for( i=0; i<njob; i++ )
1941	{
1942	myfgets( tmpname, N-1, fp );
1943	tmpseq = load1SeqWithoutName_realloc_casepreserve( fp );
1944	tmp = strlen( tmpseq );
1945	if( tmp > nlenmax ) nlenmax = tmp;
1946	if( tmp < nlenminpt ) nlenminpt = tmp;
1947	atgcnum += countATGC( tmpseq, &nsite );
1948	total += nsite;
1949	free( tmpseq );
1950	}
1951	free( tmpname );
1952	atgcfreq = (double)atgcnum / total;
1953	// fprintf( stderr, "##### atgcfreq = %f\n", atgcfreq );
1954	if( dorp == NOTSPECIFIED )
1955	{
1956	if( atgcfreq > 0.75 )
1957	{
1958	dorp = 'd';
1959	upperCase = -1;
1960	}
1961	else
1962	{
1963	dorp = 'p';
1964	upperCase = 0;
1965	}
1966	}
1967	}
1968
1969	void getnumlen_nogap( FILE fp, int nlenminpt )
1970	{
1971	int total;
1972	int nsite = 0;
1973	int atgcnum;
1974	int i, tmp;
1975	char tmpseq, tmpname;
1976	double atgcfreq;
1977	tmpname = AllocateCharVec( N );
1978	njob = countKUorWA( fp );
1979	searchKUorWA( fp );
1980	nlenmax = 0;
1981	*nlenminpt = 99999999;
1982	atgcnum = 0;
1983	total = 0;
1984	for( i=0; i<njob; i++ )
1985	{
1986	myfgets( tmpname, N-1, fp );
1987	tmpseq = load1SeqWithoutName_realloc( fp );
1988	tmp = countnogaplen( tmpseq );
1989	if( tmp > nlenmax ) nlenmax = tmp;
1990	if( tmp < nlenminpt ) nlenminpt = tmp;
1991	atgcnum += countATGC( tmpseq, &nsite );
1992	total += nsite;
1993	free( tmpseq );
1994	}
1995	free( tmpname );
1996	atgcfreq = (double)atgcnum / total;
1997	fprintf( stderr, "##### atgcfreq = %f\n", atgcfreq );
1998	if( dorp == NOTSPECIFIED )
1999	{
2000	if( atgcfreq > 0.75 )
2001	{
2002	dorp = 'd';
2003	upperCase = -1;
2004	}
2005	else
2006	{
2007	dorp = 'p';
2008	upperCase = 0;
2009	}
2010	}
2011	}
2012
2013
2014	void getnumlen_nogap_outallreg( FILE fp, int nlenminpt )
2015	{
2016	int total;
2017	int nsite = 0;
2018	int atgcnum;
2019	int i, tmp;
2020	char tmpseq, tmpname;
2021	double atgcfreq;
2022	tmpname = AllocateCharVec( N );
2023	njob = countKUorWA( fp );
2024	searchKUorWA( fp );
2025	nlenmax = 0;
2026	*nlenminpt = 99999999;
2027	atgcnum = 0;
2028	total = 0;
2029	for( i=0; i<njob; i++ )
2030	{
2031	myfgets( tmpname, N-1, fp );
2032	fprintf( stdout, "%s\n", tmpname );
2033	tmpseq = load1SeqWithoutName_realloc_casepreserve( fp );
2034	tmp = countnogaplen( tmpseq );
2035	fprintf( stdout, "%d\n", tmp );
2036	if( tmp > nlenmax ) nlenmax = tmp;
2037	if( tmp < nlenminpt ) nlenminpt = tmp;
2038	atgcnum += countATGC( tmpseq, &nsite );
2039	total += nsite;
2040	free( tmpseq );
2041	}
2042	free( tmpname );
2043	atgcfreq = (double)atgcnum / total;
2044	fprintf( stderr, "##### atgcfreq = %f\n", atgcfreq );
2045	if( dorp == NOTSPECIFIED )
2046	{
2047	if( atgcfreq > 0.75 )
2048	{
2049	dorp = 'd';
2050	upperCase = -1;
2051	}
2052	else
2053	{
2054	dorp = 'p';
2055	upperCase = 0;
2056	}
2057	}
2058	}
2059
2060	static void escapehtml( char res, char ori, int maxlen )
2061	{
2062	char *res0 = res;
2063	while( *ori )
2064	{
2065	if( *ori == '<' )
2066	{
2067	strcpy( res, "<" );
2068	res += 3;
2069	}
2070	else if( *ori == '>' )
2071	{
2072	strcpy( res, ">" );
2073	res += 3;
2074	}
2075	else if( *ori == '&' )
2076	{
2077	strcpy( res, "&" );
2078	res += 4;
2079	}
2080	else if( *ori == '"' )
2081	{
2082	strcpy( res, """ );
2083	res += 5;
2084	}
2085	else if( *ori == ' ' )
2086	{
2087	strcpy( res, " " );
2088	res += 5;
2089	}
2090	else
2091	{
2092	res = ori;
2093	}
2094	res++;
2095	ori++;
2096
2097	if( res - res0 -10 > N ) break;
2098	}
2099	*res = 0;
2100	}
2101
2102	void getnumlen_nogap_outallreg_web( FILE fp, FILE ofp, int nlenminpt, int isalignedpt )
2103	{
2104	int total;
2105	int nsite = 0;
2106	int atgcnum;
2107	int alnlen = 0, alnlen_prev;
2108	int i, tmp, lennormalchar;
2109	char tmpseq, tmpname, *tmpname2;
2110	double atgcfreq;
2111	tmpname = AllocateCharVec( N );
2112	tmpname2 = AllocateCharVec( N );
2113	njob = countKUorWA( fp );
2114	searchKUorWA( fp );
2115	nlenmax = 0;
2116	*nlenminpt = 99999999;
2117	atgcnum = 0;
2118	total = 0;
2119	alnlen_prev = -1;
2120	*isalignedpt = 1;
2121	for( i=0; i<njob; i++ )
2122	{
2123	myfgets( tmpname, N-1, fp );
2124	tmpname2[0] = tmpname[0];
2125	escapehtml( tmpname2+1, tmpname+1, N );
2126	// fprintf( stdout, "%s\n", tmpname );
2127	// fprintf( stdout, "%s\n", tmpname2 );
2128	// exit(1);
2129	tmpseq = load1SeqWithoutName_realloc_casepreserve( fp );
2130	tmp = countnogaplen( tmpseq );
2131	// fprintf( stdout, "%d\n", tmp );
2132	if( tmp > nlenmax ) nlenmax = tmp;
2133	if( tmp < nlenminpt ) nlenminpt = tmp;
2134	atgcnum += countATGC( tmpseq, &nsite );
2135	total += nsite;
2136
2137	alnlen = strlen( tmpseq );
2138	// fprintf( stdout, "##### alnlen, alnlen_prev = %d, %d\n", alnlen, alnlen_prev );
2139	if( i>0 && alnlen_prev != alnlen ) *isalignedpt = 0;
2140	alnlen_prev = alnlen;
2141
2142	atgcfreq = (double)atgcnum / total;
2143	// fprintf( stderr, "##### atgcfreq = %f\n", atgcfreq );
2144	// if( dorp == NOTSPECIFIED ) // you kentou
2145	{
2146	if( atgcfreq > 0.75 )
2147	{
2148	dorp = 'd';
2149	upperCase = -1;
2150	}
2151	else
2152	{
2153	dorp = 'p';
2154	upperCase = 0;
2155	}
2156	}
2157
2158	if( dorp == 'd' ) lennormalchar = countnormalletters( tmpseq, "atgcuATGCU" );
2159	else lennormalchar = countnormalletters( tmpseq, "ARNDCQEGHILKMFPSTWYVarndcqeghilkmfpstwyv" );
2160	free( tmpseq );
2161
2162	fprintf( ofp, " <label for='s%d'><span id='ss%d'><input type='checkbox' id='s%d' name='s%d' checked></span> <input type='text' class='ll' id='ll%d' style='display:none' size='6' value='%d' readonly='readonly'>%s</label>\n", i, i, i, i, i, lennormalchar, tmpname2 );
2163	fprintf( ofp, "<span id='t%d-0' style='display:none'>", i );
2164	fprintf( ofp, " <a href='javascript:void(0)' onclick='ddcycle(this.form,\"t%d\")'>+reg</a>", i );
2165	fprintf( ofp, " Begin:<input type='text' name='b%d-0' size='8' value='1' class='ie'> End:<input type='text' name='e%d-0' size='8' value='%d' class='ie'>", i, i, tmp );
2166	if( dorp == 'd' ) fprintf( ofp, " <label for='r%d-0'><input type='checkbox' name='r%d-0' id='r%d-0'>Reverse</label>", i, i, i );
2167	// fprintf( ofp, " Sequence Length:<input type='text' name='l%d' size='8' value='%d' readonly='readonly'>", i, tmp );
2168	fprintf( ofp, "\n</span>" );
2169	fprintf( ofp, "<span id='t%d-1' style='display:none'>", i );
2170	fprintf( ofp, " Begin:<input type='text' name='b%d-1' size='8' value='' class='ie'> End:<input type='text' name='e%d-1' size='8' value='' class='ie'>", i, i );
2171	if( dorp == 'd' ) fprintf( ofp, " <label for='r%d-1'><input type='checkbox' name='r%d-1' id='r%d-1'>Reverse</label>", i, i, i );
2172	fprintf( ofp, "\n</span>" );
2173	fprintf( ofp, "<span id='t%d-2' style='display:none'>", i );
2174	fprintf( ofp, " Begin:<input type='text' name='b%d-2' size='8' value='' class='ie'> End:<input type='text' name='e%d-2' size='8' value='' class='ie'>", i, i );
2175	if( dorp == 'd' ) fprintf( ofp, " <label for='r%d-2'><input type='checkbox' name='r%d-2' id='r%d-2'>Reverse</label>", i, i, i );
2176	fprintf( ofp, "\n</span>" );
2177	fprintf( ofp, "<span id='t%d-3' style='display:none'>", i );
2178	fprintf( ofp, " Begin:<input type='text' name='b%d-3' size='8' value='' class='ie'> End:<input type='text' name='e%d-3' size='8' value='' class='ie'>", i, i );
2179	if( dorp == 'd' ) fprintf( ofp, " <label for='r%d-3'><input type='checkbox' name='r%d-3' id='r%d-3'>Reverse</label>", i, i, i );
2180	fprintf( ofp, "\n</span>" );
2181	fprintf( ofp, "<span id='t%d-4' style='display:none'>", i );
2182	fprintf( ofp, " Begin:<input type='text' name='b%d-4' size='8' value='' class='ie'> End:<input type='text' name='e%d-4' size='8' value='' class='ie'>", i, i );
2183	if( dorp == 'd' ) fprintf( ofp, " <label for='r%d-4'><input type='checkbox' name='r%d-4' id='r%d-4'>Reverse</label>", i, i, i );
2184	fprintf( ofp, "\n</span>" );
2185	}
2186	free( tmpname );
2187	free( tmpname2 );
2188	atgcfreq = (double)atgcnum / total;
2189	fprintf( stderr, "##### atgcfreq = %f\n", atgcfreq );
2190	// if( dorp == NOTSPECIFIED ) // you kentou
2191	{
2192	if( atgcfreq > 0.75 )
2193	{
2194	dorp = 'd';
2195	upperCase = -1;
2196	}
2197	else
2198	{
2199	dorp = 'p';
2200	upperCase = 0;
2201	}
2202	}
2203	fprintf( ofp, "\n" );
2204	if( *isalignedpt )
2205	{
2206	fprintf( ofp, "<span id='tall-0' style='display:none'>" );
2207	fprintf( ofp, "Cut the alignment\n" );
2208	fprintf( ofp, " <a href='javascript:void(0)' onclick='ddcycle(this.form,\"tall\")'>+reg</a>" );
2209	fprintf( ofp, " Begin:<input type='text' name='ball-0' size='8' value='1'> End:<input type='text' name='eall-0' size='8' value='%d'>", alnlen );
2210	if( dorp == 'd' ) fprintf( ofp, " <label for='rall-0'><input type='checkbox' name='rall-0' id='rall-0'>Reverse</label>" );
2211	fprintf( ofp, " Alignment length:<input type='text' name='lall' size='8' value='%d' readonly='readonly'>", alnlen );
2212	fprintf( ofp, "\n</span>" );
2213	fprintf( ofp, "<span id='tall-1' style='display:none'>" );
2214	fprintf( ofp, " Begin:<input type='text' name='ball-1' size='8' value=''> End:<input type='text' name='eall-1' size='8' value=''>" );
2215	if( dorp == 'd' ) fprintf( ofp, " <label for='rall-1'><input type='checkbox' name='rall-1' id='rall-1'>Reverse</label>" );
2216	fprintf( ofp, "\n</span>" );
2217	fprintf( ofp, "<span id='tall-2' style='display:none'>" );
2218	fprintf( ofp, " Begin:<input type='text' name='ball-2' size='8' value=''> End:<input type='text' name='eall-1' size='8' value=''>" );
2219	if( dorp == 'd' ) fprintf( ofp, " <label for='rall-2'><input type='checkbox' name='rall-2' id='rall-2'>Reverse</label>" );
2220	fprintf( ofp, "\n</span>" );
2221	fprintf( ofp, "<span id='tall-3' style='display:none'>" );
2222	fprintf( ofp, " Begin:<input type='text' name='ball-3' size='8' value=''> End:<input type='text' name='eall-1' size='8' value=''>" );
2223	if( dorp == 'd' ) fprintf( ofp, " <label for='rall-3'><input type='checkbox' name='rall-3' id='rall-3'>Reverse</label>" );
2224	fprintf( ofp, "\n</span>" );
2225	fprintf( ofp, "<span id='tall-4' style='display:none'>" );
2226	fprintf( ofp, " Begin:<input type='text' name='ball-4' size='8' value=''> End:<input type='text' name='eall-1' size='8' value=''>" );
2227	if( dorp == 'd' ) fprintf( ofp, " <label for='rall-4'><input type='checkbox' name='rall-4' id='rall-4'>Reverse</label>" );
2228	fprintf( ofp, "\n</span>" );
2229	}
2230
2231	}
2232
2233	void getnumlen( FILE *fp )
2234	{
2235	int total;
2236	int nsite = 0;
2237	int atgcnum;
2238	int i, tmp;
2239	char *tmpseq;
2240	char *tmpname;
2241	double atgcfreq;
2242	tmpname = AllocateCharVec( N );
2243	njob = countKUorWA( fp );
2244	searchKUorWA( fp );
2245	nlenmax = 0;
2246	atgcnum = 0;
2247	total = 0;
2248	for( i=0; i<njob; i++ )
2249	{
2250	myfgets( tmpname, N-1, fp );
2251	tmpseq = load1SeqWithoutName_realloc( fp );
2252	tmp = strlen( tmpseq );
2253	if( tmp > nlenmax ) nlenmax = tmp;
2254	atgcnum += countATGC( tmpseq, &nsite );
2255	total += nsite;
2256	// fprintf( stderr, "##### total = %d\n", total );
2257	free( tmpseq );
2258	}
2259
2260
2261	atgcfreq = (double)atgcnum / total;
2262	// fprintf( stderr, "##### atgcfreq = %f\n", atgcfreq );
2263	if( dorp == NOTSPECIFIED )
2264	{
2265	if( atgcfreq > 0.75 )
2266	{
2267	dorp = 'd';
2268	upperCase = -1;
2269	}
2270	else
2271	{
2272	dorp = 'p';
2273	upperCase = 0;
2274	}
2275	}
2276	free( tmpname );
2277	}
2278
2279
2280
2281	void WriteGapFill( FILE fp, int locnjob, char name[][B], int nlen[M], char *aseq )
2282	{
2283	static char b[N];
2284	int i, j;
2285	int nalen[M];
2286	static char gap[N];
2287	static char buff[N];
2288
2289	#if IODEBUG
2290	fprintf( stderr, "IMAKARA KAKU\n" );
2291	#endif
2292	nlenmax = 0;
2293	for( i=0; i<locnjob; i++ )
2294	{
2295	int len = strlen( aseq[i] );
2296	if( nlenmax < len ) nlenmax = len;
2297	}
2298
2299	for( i=0; i<nlenmax; i++ ) gap[i] = '-';
2300	gap[nlenmax] = 0;
2301
2302	fprintf( fp, "%5d", locnjob );
2303	fprintf( fp, "\n" );
2304
2305	for( i=0; i<locnjob; i++ )
2306	{
2307	strcpy( buff, aseq[i] );
2308	strncat( buff, gap, nlenmax-strlen( aseq[i] ) );
2309	buff[nlenmax] = 0;
2310	nalen[i] = strlen( buff );
2311	fprintf( fp, "%s\n", name[i] );
2312	fprintf( fp, "%5d\n", nalen[i] );
2313	for( j=0; j<nalen[i]; j=j+C )
2314	{
2315	strncpy_caseC( b, buff+j, C ); b[C] = 0;
2316	fprintf( fp, "%s\n",b );
2317	}
2318	}
2319	#if DEBUG
2320	fprintf( stderr, "nalen[0] = %d\n", nalen[0] );
2321	#endif
2322	#if IODEBUG
2323	fprintf( stderr, "KAKIOWATTA\n" );
2324	#endif
2325	}
2326
2327	void writeDataforgaln( FILE fp, int locnjob, char name, int nlen, char **aseq )
2328	{
2329	int i, j;
2330	int nalen;
2331
2332	for( i=0; i<locnjob; i++ )
2333	{
2334	nalen = strlen( aseq[i] );
2335	fprintf( fp, ">%s\n", name[i]+1 );
2336	for( j=0; j<nalen; j=j+C )
2337	{
2338	#if 0
2339	strncpy( b, aseq[i]+j, C ); b[C] = 0;
2340	fprintf( fp, "%s\n",b );
2341	#else
2342	fprintf( fp, "%.*s\n", C, aseq[i]+j );
2343	#endif
2344	}
2345	}
2346	}
2347
2348	void writeData_pointer( FILE fp, int locnjob, char name, int nlen, char **aseq )
2349	{
2350	int i, j;
2351	int nalen;
2352
2353	for( i=0; i<locnjob; i++ )
2354	{
2355	#if DEBUG
2356	fprintf( stderr, "i = %d in writeData\n", i );
2357	#endif
2358	nalen = strlen( aseq[i] );
2359	fprintf( fp, ">%s\n", name[i]+1 );
2360	for( j=0; j<nalen; j=j+C )
2361	{
2362	#if 0
2363	strncpy( b, aseq[i]+j, C ); b[C] = 0;
2364	fprintf( fp, "%s\n",b );
2365	#else
2366	fprintf( fp, "%.*s\n", C, aseq[i]+j );
2367	#endif
2368	}
2369	}
2370	}
2371
2372	void writeData( FILE fp, int locnjob, char name[][B], int nlen[], char *aseq )
2373	{
2374	int i, j;
2375	int nalen;
2376
2377	for( i=0; i<locnjob; i++ )
2378	{
2379	#if DEBUG
2380	fprintf( stderr, "i = %d in writeData\n", i );
2381	#endif
2382	nalen = strlen( aseq[i] );
2383	fprintf( fp, ">%s\n", name[i]+1 );
2384	for( j=0; j<nalen; j=j+C )
2385	{
2386	#if 0
2387	strncpy( b, aseq[i]+j, C ); b[C] = 0;
2388	fprintf( fp, "%s\n",b );
2389	#else
2390	fprintf( fp, "%.*s\n", C, aseq[i]+j );
2391	#endif
2392	}
2393	}
2394	}
2395
2396
2397	void write1seq( FILE fp, char aseq )
2398	{
2399	int j;
2400	int nalen;
2401
2402	nalen = strlen( aseq );
2403	for( j=0; j<nalen; j=j+C )
2404	fprintf( fp, "%.*s\n", C, aseq+j );
2405	}
2406
2407	void readhat2_floathalf_part_pointer( FILE fp, int nseq, int nadd, char name, float *mtx )
2408	{
2409	int i, j, nseq0, norg;
2410	char b[B];
2411
2412	fgets( b, B, fp );
2413	fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 )
2414	{
2415	fprintf( stderr, "%d != %d\n", nseq, nseq0 );
2416	ErrorExit( "hat2 is wrong." );
2417	}
2418	fgets( b, B, fp );
2419	for( i=0; i<nseq; i++ )
2420	{
2421	#if 0
2422	getaline_fp_eof( b, B, fp );
2423	#else
2424	myfgets( b, B-2, fp );
2425	#endif
2426	#if 0
2427	j = MIN( strlen( b+6 ), 10 );
2428	if( strncmp( name[i], b+6 , j ) )
2429	{
2430	fprintf( stderr, "Error in hat2\n" );
2431	fprintf( stderr, "%s != %s\n", b, name[i] );
2432	exit( 1 );
2433	}
2434	#endif
2435	}
2436	norg = nseq-nadd;
2437	for( i=0; i<norg; i++ ) for( j=0; j<nadd; j++ )
2438	{
2439	mtx[i][j] = ( input_new( fp, D ) );
2440	}
2441	}
2442
2443	void readhat2_floathalf_pointer( FILE fp, int nseq, char name, float *mtx )
2444	{
2445	int i, j, nseq0;
2446	char b[B];
2447
2448	fgets( b, B, fp );
2449	fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 )
2450	{
2451	fprintf( stderr, "%d != %d\n", nseq, nseq0 );
2452	ErrorExit( "hat2 is wrong." );
2453	}
2454	fgets( b, B, fp );
2455	for( i=0; i<nseq; i++ )
2456	{
2457	#if 0
2458	getaline_fp_eof( b, B, fp );
2459	#else
2460	myfgets( b, B-2, fp );
2461	#endif
2462	#if 0
2463	j = MIN( strlen( b+6 ), 10 );
2464	if( strncmp( name[i], b+6 , j ) )
2465	{
2466	fprintf( stderr, "Error in hat2\n" );
2467	fprintf( stderr, "%s != %s\n", b, name[i] );
2468	exit( 1 );
2469	}
2470	#endif
2471	}
2472	for( i=0; i<nseq-1; i++ ) for( j=i+1; j<nseq; j++ )
2473	{
2474	mtx[i][j-i] = ( input_new( fp, D ) );
2475	}
2476	}
2477	void readhat2_floathalf( FILE fp, int nseq, char name[M][B], float *mtx )
2478	{
2479	int i, j, nseq0;
2480	char b[B];
2481
2482	fgets( b, B, fp );
2483	fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 ) ErrorExit( "hat2 is wrong." );
2484	fgets( b, B, fp );
2485	for( i=0; i<nseq; i++ )
2486	{
2487	#if 0
2488	getaline_fp_eof( b, B, fp );
2489	#else
2490	myfgets( b, B-2, fp );
2491	#endif
2492	#if 0
2493	j = MIN( strlen( b+6 ), 10 );
2494	if( strncmp( name[i], b+6 , j ) )
2495	{
2496	fprintf( stderr, "Error in hat2\n" );
2497	fprintf( stderr, "%s != %s\n", b, name[i] );
2498	exit( 1 );
2499	}
2500	#endif
2501	}
2502	for( i=0; i<nseq-1; i++ ) for( j=i+1; j<nseq; j++ )
2503	{
2504	mtx[i][j-i] = ( input_new( fp, D ) );
2505	}
2506	}
2507	void readhat2_float( FILE fp, int nseq, char name[M][B], float *mtx )
2508	{
2509	int i, j, nseq0;
2510	char b[B];
2511
2512	fgets( b, B, fp );
2513	fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 ) ErrorExit( "hat2 is wrong." );
2514	fgets( b, B, fp );
2515	for( i=0; i<nseq; i++ )
2516	{
2517	#if 0
2518	getaline_fp_eof( b, B, fp );
2519	#else
2520	myfgets( b, B-2, fp );
2521	#endif
2522	#if 0
2523	j = MIN( strlen( b+6 ), 10 );
2524	if( strncmp( name[i], b+6 , j ) )
2525	{
2526	fprintf( stderr, "Error in hat2\n" );
2527	fprintf( stderr, "%s != %s\n", b, name[i] );
2528	exit( 1 );
2529	}
2530	#endif
2531	}
2532	for( i=0; i<nseq-1; i++ ) for( j=i+1; j<nseq; j++ )
2533	{
2534	mtx[i][j] = ( input_new( fp, D ) );
2535	}
2536	}
2537	void readhat2_int( FILE fp, int nseq, char name[M][B], int *mtx )
2538	{
2539	int i, j, nseq0;
2540	char b[B];
2541
2542	fgets( b, B, fp );
2543	fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 ) ErrorExit( "hat2 is wrong." );
2544	fgets( b, B, fp );
2545	for( i=0; i<nseq; i++ )
2546	{
2547	#if 0
2548	getaline_fp_eof( b, B, fp );
2549	#else
2550	myfgets( b, B-2, fp );
2551	#endif
2552	#if 0
2553	j = MIN( strlen( b+6 ), 10 );
2554	if( strncmp( name[i], b+6 , j ) )
2555	{
2556	fprintf( stderr, "Error in hat2\n" );
2557	fprintf( stderr, "%s != %s\n", b, name[i] );
2558	exit( 1 );
2559	}
2560	#endif
2561	}
2562	for( i=0; i<nseq-1; i++ ) for( j=i+1; j<nseq; j++ )
2563	{
2564	mtx[i][j] = (int)( input_new( fp, D ) * INTMTXSCALE + 0.5 );
2565	}
2566	}
2567
2568	void readhat2_pointer( FILE fp, int nseq, char name, double *mtx )
2569	{
2570	int i, j, nseq0;
2571	char b[B];
2572
2573	fgets( b, B, fp );
2574	fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 ) ErrorExit( "hat2 is wrong." );
2575	fgets( b, B, fp );
2576	for( i=0; i<nseq; i++ )
2577	{
2578	#if 0
2579	getaline_fp_eof( b, B, fp );
2580	#else
2581	myfgets( b, B-2, fp );
2582	#endif
2583	#if 0
2584	j = MIN( strlen( b+6 ), 10 );
2585	if( strncmp( name[i], b+6 , j ) )
2586	{
2587	fprintf( stderr, "Error in hat2\n" );
2588	fprintf( stderr, "%s != %s\n", b, name[i] );
2589	exit( 1 );
2590	}
2591	#endif
2592	}
2593	for( i=0; i<nseq-1; i++ ) for( j=i+1; j<nseq; j++ )
2594	{
2595	mtx[i][j] = (double)input_new( fp, D);
2596	}
2597	}
2598	void readhat2( FILE fp, int nseq, char name[M][B], double *mtx )
2599	{
2600	int i, j, nseq0;
2601	char b[B];
2602
2603	fgets( b, B, fp );
2604	fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 ) ErrorExit( "hat2 is wrong." );
2605	fgets( b, B, fp );
2606	for( i=0; i<nseq; i++ )
2607	{
2608	#if 0
2609	getaline_fp_eof( b, B, fp );
2610	#else
2611	myfgets( b, B-2, fp );
2612	#endif
2613	#if 0
2614	j = MIN( strlen( b+6 ), 10 );
2615	if( strncmp( name[i], b+6 , j ) )
2616	{
2617	fprintf( stderr, "Error in hat2\n" );
2618	fprintf( stderr, "%s != %s\n", b, name[i] );
2619	exit( 1 );
2620	}
2621	#endif
2622	}
2623	for( i=0; i<nseq-1; i++ ) for( j=i+1; j<nseq; j++ )
2624	{
2625	mtx[i][j] = (double)input_new( fp, D);
2626	}
2627	}
2628
2629	void WriteFloatHat2_pointer_halfmtx( FILE hat2p, int locnjob, char name, float *mtx )
2630	{
2631	int i, j, ijsa;
2632	double max = 0.0;
2633	for( i=0; i<locnjob-1; i++ ) for( j=1; j<locnjob-i; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
2634
2635	fprintf( hat2p, "%5d\n", 1 );
2636	fprintf( hat2p, "%5d\n", locnjob );
2637	fprintf( hat2p, " %#6.3f\n", max * 2.5 );
2638
2639	for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
2640	for( i=0; i<locnjob; i++ )
2641	{
2642	for( j=i+1; j<njob; j++ )
2643	{
2644	fprintf( hat2p, "%#6.3f", mtx[i][j-i] );
2645	ijsa = j-i;
2646	if( ijsa % 12 == 0 \|\| ijsa == locnjob-i-1 ) fprintf( hat2p, "\n" );
2647	}
2648	}
2649	}
2650
2651	void WriteFloatHat2_pointer( FILE hat2p, int locnjob, char name, float *mtx )
2652	{
2653	int i, j;
2654	double max = 0.0;
2655	for( i=0; i<locnjob-1; i++ ) for( j=1; j<locnjob-i; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
2656
2657	fprintf( hat2p, "%5d\n", 1 );
2658	fprintf( hat2p, "%5d\n", locnjob );
2659	fprintf( hat2p, " %#6.3f\n", max * 2.5 );
2660
2661	for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
2662	for( i=0; i<locnjob; i++ )
2663	{
2664	for( j=1; j<locnjob-i; j++ )
2665	{
2666	fprintf( hat2p, "%#6.3f", mtx[i][j] );
2667	if( j % 12 == 0 \|\| j == locnjob-i-1 ) fprintf( hat2p, "\n" );
2668	}
2669	}
2670	}
2671
2672	void WriteFloatHat2( FILE hat2p, int locnjob, char name[M][B], float *mtx )
2673	{
2674	int i, j;
2675	double max = 0.0;
2676	for( i=0; i<locnjob-1; i++ ) for( j=1; j<locnjob-i; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
2677
2678	fprintf( hat2p, "%5d\n", 1 );
2679	fprintf( hat2p, "%5d\n", locnjob );
2680	fprintf( hat2p, " %#6.3f\n", max * 2.5 );
2681
2682	for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
2683	for( i=0; i<locnjob; i++ )
2684	{
2685	for( j=1; j<locnjob-i; j++ )
2686	{
2687	fprintf( hat2p, "%#6.3f", mtx[i][j] );
2688	if( j % 12 == 0 \|\| j == locnjob-i-1 ) fprintf( hat2p, "\n" );
2689	}
2690	}
2691	}
2692
2693	void WriteHat2_int( FILE hat2p, int locnjob, char name[M][B], int *mtx )
2694	{
2695	int i, j;
2696	double max = 0.0;
2697	for( i=0; i<locnjob-1; i++ ) for( j=i+1; j<locnjob; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
2698	max /= INTMTXSCALE;
2699
2700	fprintf( hat2p, "%5d\n", 1 );
2701	fprintf( hat2p, "%5d\n", locnjob );
2702	fprintf( hat2p, " %#6.3f\n", max * 2.5 );
2703
2704	for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
2705	for( i=0; i<locnjob-1; i++ )
2706	{
2707	for( j=i+1; j<locnjob; j++ )
2708	{
2709	fprintf( hat2p, "%#6.3f", (float)mtx[i][j] / INTMTXSCALE );
2710	if( (j-i) % 12 == 0 \|\| j == locnjob-1 ) fprintf( hat2p, "\n" );
2711	}
2712	}
2713	}
2714
2715	void WriteHat2_part_pointer( FILE hat2p, int locnjob, int nadd, char name, double *mtx )
2716	{
2717	int i, j;
2718	int norg = locnjob-nadd;
2719	double max = 0.0;
2720	// for( i=0; i<locnjob-1; i++ ) for( j=i+1; j<locnjob; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
2721
2722	fprintf( hat2p, "%5d\n", 1 );
2723	fprintf( hat2p, "%5d\n", locnjob );
2724	fprintf( hat2p, " %#6.3f\n", max * 2.5 );
2725
2726	for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
2727	for( i=0; i<norg; i++ )
2728	{
2729	for( j=0; j<nadd; j++ )
2730	{
2731	fprintf( hat2p, "%#6.3f", mtx[i][j] );
2732	if( (j+1) % 12 == 0 \|\| j == nadd-1 ) fprintf( hat2p, "\n" );
2733	}
2734	}
2735	}
2736
2737	void WriteHat2_pointer( FILE hat2p, int locnjob, char name, double *mtx )
2738	{
2739	int i, j;
2740	double max = 0.0;
2741	for( i=0; i<locnjob-1; i++ ) for( j=i+1; j<locnjob; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
2742
2743	fprintf( hat2p, "%5d\n", 1 );
2744	fprintf( hat2p, "%5d\n", locnjob );
2745	fprintf( hat2p, " %#6.3f\n", max * 2.5 );
2746
2747	for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
2748	for( i=0; i<locnjob-1; i++ )
2749	{
2750	for( j=i+1; j<locnjob; j++ )
2751	{
2752	fprintf( hat2p, "%#6.3f", mtx[i][j] );
2753	if( (j-i) % 12 == 0 \|\| j == locnjob-1 ) fprintf( hat2p, "\n" );
2754	}
2755	}
2756	}
2757
2758	void WriteHat2( FILE hat2p, int locnjob, char name[M][B], double *mtx )
2759	{
2760	int i, j;
2761	double max = 0.0;
2762	for( i=0; i<locnjob-1; i++ ) for( j=i+1; j<locnjob; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
2763
2764	fprintf( hat2p, "%5d\n", 1 );
2765	fprintf( hat2p, "%5d\n", locnjob );
2766	fprintf( hat2p, " %#6.3f\n", max * 2.5 );
2767
2768	for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
2769	for( i=0; i<locnjob-1; i++ )
2770	{
2771	for( j=i+1; j<locnjob; j++ )
2772	{
2773	fprintf( hat2p, "%#6.3f", mtx[i][j] );
2774	if( (j-i) % 12 == 0 \|\| j == locnjob-1 ) fprintf( hat2p, "\n" );
2775	}
2776	}
2777	}
2778
2779	#if 0
2780	void WriteHat2plain( FILE hat2p, int locnjob, double *mtx )
2781	{
2782	int i, j, ilim;
2783
2784	ilim = locnjob-1;
2785	for( i=0; i<ilim; i++ )
2786	{
2787	fprintf( hat2p, "%d-%d d=%.3f\n", i+1, i+1, 0.0 );
2788	for( j=i+1; j<locnjob; j++ )
2789	{
2790	fprintf( hat2p, "%d-%d d=%.3f\n", i+1, j+1, mtx[i][j] );
2791	}
2792	}
2793	}
2794	#endif
2795
2796	int ReadFasta_sub( FILE fp, double dis, int nseq, char name[M][B] )
2797	{
2798	int i, count=0;
2799	char b[B];
2800	int junban[M];
2801
2802	count = 0;
2803	for( i=0; i<10000000 && count<nseq; i++ )
2804	{
2805	fgets( b, B-1, fp );
2806	if( !strncmp( "+==========+", b, 12 ) )
2807	{
2808	junban[count] = atoi( b+12 );
2809	count++;
2810	}
2811	}
2812
2813	for( i=0; i<nseq; i++ ) dis[i] = 0.0;
2814	count = 0;
2815	for( i=0; i<100000 && count<nseq; i++ )
2816	{
2817	if( fgets( b, B-1, fp ) ) break;
2818	if( !strncmp( name[junban[count]], b, 20 ) )
2819	{
2820	fgets( b, B-1, fp );
2821	dis[junban[count]] = atof( b );
2822	count++;
2823	}
2824	}
2825	return 0;
2826	}
2827
2828
2829	int ReadSsearch( FILE fp, double dis, int nseq, char name[M][B] )
2830	{
2831	int i, count=0;
2832	char b[B];
2833	int junban[M];
2834	int opt;
2835
2836	count = 0;
2837	for( i=0; i<10000000 && count<nseq; i++ )
2838	{
2839	fgets( b, B-1, fp );
2840	if( !strncmp( "+==========+", b, 12 ) )
2841	{
2842	junban[count] = atoi( b+12 );
2843	sscanf( b+75, "%d", &opt );
2844	dis[junban[count]] = (double)opt;
2845	count++;
2846	}
2847	}
2848
2849	/*
2850	count = 0;
2851	for( i=0; i<100000 && count<nseq; i++ )
2852	{
2853	fgets( b, B-1, fp );
2854	if( !strncmp( name[junban[count]], b, 20 ) )
2855	{
2856	dis[junban[count]] = atof( b+65 );
2857	count++;
2858	}
2859	}
2860	*/
2861	return 0;
2862	}
2863
2864	int ReadBlastm7_avscore( FILE fp, double dis, int nin )
2865	{
2866	int count=0;
2867	char b[B];
2868	char *pt;
2869	int *junban;
2870	double score, sumscore;
2871	double len, sumlen;
2872	int qstart, qend, tstart, tend;
2873	double scorepersite;
2874	static char qal[N], tal[N], al[N];
2875	int nlocalhom;
2876
2877	junban = calloc( nin, sizeof( int ) );
2878
2879	count = 0;
2880	sumscore = 0.0;
2881	sumlen = 0.0;
2882	score = 0.0;
2883	len = 0.0;
2884	scorepersite = 0.0; // by D.Mathog, a guess
2885	while( 1 )
2886	{
2887
2888	if( feof( fp ) ) break;
2889
2890	while( fgets( b, B-1, fp ) )
2891	{
2892	if( !strncmp( " <Hit_def>", b, 19 ) \|\| !strncmp( " <Hsp_num>", b, 23 ) ) break;
2893	}
2894
2895	if( !strncmp( " <Hit_def>", b, 19 ) )
2896	{
2897	junban[count] = atoi( b+31 );
2898	nlocalhom = 0;
2899	}
2900
2901
2902	while( fgets( b, B-1, fp ) )
2903	if( !strncmp( " <Hsp_score>", b, 25 ) ) break;
2904	pt = b + 25;
2905	score = atof( pt );
2906	sumscore += score;
2907
2908
2909	while( fgets( b, B-1, fp ) )
2910	if( !strncmp( " <Hsp_query-from>", b, 30 ) ) break;
2911	pt = b + 30;
2912	qstart = atoi( pt ) - 1;
2913
2914
2915	while( fgets( b, B-1, fp ) )
2916	if( !strncmp( " <Hsp_query-to>", b, 28 ) ) break;
2917	pt = b + 28;
2918	qend = atoi( pt ) - 1;
2919
2920
2921	while( fgets( b, B-1, fp ) )
2922	if( !strncmp( " <Hsp_hit-from>", b, 28 ) ) break;
2923	pt = b + 28;
2924	tstart = atoi( pt ) - 1;
2925
2926
2927	while( fgets( b, B-1, fp ) )
2928	if( !strncmp( " <Hsp_hit-to>", b, 26 ) ) break;
2929	pt = b + 26;
2930	tend = atoi( pt ) - 1;
2931
2932
2933	while( fgets( b, B-1, fp ) )
2934	if( !strncmp( " <Hsp_align-len>", b, 29 ) ) break;
2935	pt = b + 29;
2936	len = atoi( pt );
2937	sumlen += len;
2938
2939
2940	while( fgets( al, N-100, fp ) )
2941	if( !strncmp( " <Hsp_qseq>", al, 24 ) ) break;
2942
2943	strcpy( qal, al+24 );
2944	pt = qal;
2945	while( *++pt != '<' )
2946	;
2947	*pt = 0;
2948
2949
2950	while( fgets( al, N-100, fp ) )
2951	if( !strncmp( " <Hsp_hseq>", al, 24 ) ) break;
2952
2953	strcpy( tal, al+24 );
2954	pt = tal;
2955	while( *++pt != '<' )
2956	;
2957	*pt = 0;
2958
2959
2960	// fprintf( stderr, "t=%d, score = %f, qstart=%d, qend=%d, tstart=%d, tend=%d, overlapaa=%d\n", junban[count], score, qstart, qend, tstart, tend, overlapaa );
2961
2962
2963	while( fgets( b, B-1, fp ) )
2964	if( !strncmp( " </Hsp>:", b, 18 ) ) break;
2965
2966
2967	fgets( b, B-1, fp );
2968
2969
2970	if( !strncmp( " </Hit_hsps>", b, 21 ) )
2971	{
2972	dis[junban[count++]] = sumscore;
2973	sumscore = 0.0;
2974	fgets( b, B-1, fp );
2975	fgets( b, B-1, fp );
2976	scorepersite = sumscore / sumlen;
2977	if( scorepersite != (int)scorepersite )
2978	{
2979	fprintf( stderr, "ERROR! sumscore=%f, sumlen=%f, and scorepersite=%f\n", sumscore, sumlen, scorepersite );
2980	exit( 1 );
2981	}
2982
2983	if( !strncmp( " </Iteration_hits>", b, 23 ) ) break;
2984	}
2985	}
2986
2987	free( junban );
2988
2989	return (int)scorepersite;
2990	}
2991	int ReadBlastm7_scoreonly( FILE fp, double dis, int nin )
2992	{
2993	int count=0;
2994	char b[B];
2995	char *pt;
2996	int *junban;
2997	int overlapaa;
2998	double score, sumscore;
2999	int qstart, qend, tstart, tend;
3000	static char qal[N], tal[N], al[N];
3001	int nlocalhom;
3002
3003	junban = calloc( nin, sizeof( int ) );
3004
3005	count = 0;
3006	sumscore = 0.0;
3007	score = 0.0;
3008	while( 1 )
3009	{
3010
3011	if( feof( fp ) ) break;
3012
3013	while( fgets( b, B-1, fp ) )
3014	{
3015	if( !strncmp( " <Hit_def>", b, 19 ) \|\| !strncmp( " <Hsp_num>", b, 23 ) ) break;
3016	}
3017
3018	if( !strncmp( " <Hit_def>", b, 19 ) )
3019	{
3020	junban[count] = atoi( b+31 );
3021	nlocalhom = 0;
3022	}
3023
3024
3025	while( fgets( b, B-1, fp ) )
3026	if( !strncmp( " <Hsp_score>", b, 25 ) ) break;
3027	pt = b + 25;
3028	score = atof( pt );
3029	sumscore += score;
3030
3031
3032	while( fgets( b, B-1, fp ) )
3033	if( !strncmp( " <Hsp_query-from>", b, 30 ) ) break;
3034	pt = b + 30;
3035	qstart = atoi( pt ) - 1;
3036
3037
3038	while( fgets( b, B-1, fp ) )
3039	if( !strncmp( " <Hsp_query-to>", b, 28 ) ) break;
3040	pt = b + 28;
3041	qend = atoi( pt ) - 1;
3042
3043
3044	while( fgets( b, B-1, fp ) )
3045	if( !strncmp( " <Hsp_hit-from>", b, 28 ) ) break;
3046	pt = b + 28;
3047	tstart = atoi( pt ) - 1;
3048
3049
3050	while( fgets( b, B-1, fp ) )
3051	if( !strncmp( " <Hsp_hit-to>", b, 26 ) ) break;
3052	pt = b + 26;
3053	tend = atoi( pt ) - 1;
3054
3055
3056	while( fgets( b, B-1, fp ) )
3057	if( !strncmp( " <Hsp_align-len>", b, 29 ) ) break;
3058	pt = b + 29;
3059	overlapaa = atoi( pt );
3060
3061
3062	while( fgets( al, N-100, fp ) )
3063	if( !strncmp( " <Hsp_qseq>", al, 24 ) ) break;
3064
3065	strcpy( qal, al+24 );
3066	pt = qal;
3067	while( *++pt != '<' )
3068	;
3069	*pt = 0;
3070
3071
3072	while( fgets( al, N-100, fp ) )
3073	if( !strncmp( " <Hsp_hseq>", al, 24 ) ) break;
3074
3075	strcpy( tal, al+24 );
3076	pt = tal;
3077	while( *++pt != '<' )
3078	;
3079	*pt = 0;
3080
3081
3082	// fprintf( stderr, "t=%d, score = %f, qstart=%d, qend=%d, tstart=%d, tend=%d, overlapaa=%d\n", junban[count], score, qstart, qend, tstart, tend, overlapaa );
3083
3084	// nlocalhom += addlocalhom_r( qal, tal, localhomlist+junban[count], qstart, tstart, score, overlapaa, nlocalhom );
3085
3086	while( fgets( b, B-1, fp ) )
3087	if( !strncmp( " </Hsp>:", b, 18 ) ) break;
3088
3089
3090	fgets( b, B-1, fp );
3091
3092
3093	if( !strncmp( " </Hit_hsps>", b, 21 ) )
3094	{
3095	dis[junban[count++]] = sumscore;
3096	sumscore = 0.0;
3097	fgets( b, B-1, fp );
3098	fgets( b, B-1, fp );
3099	if( !strncmp( " </Iteration_hits>", b, 23 ) ) break;
3100	}
3101	}
3102
3103	free( junban );
3104
3105	return count;
3106	}
3107
3108	int ReadBlastm7( FILE fp, double dis, int qmem, char *name, LocalHom localhomlist )
3109	{
3110	int count=0;
3111	char b[B];
3112	char *pt;
3113	static int junban[M];
3114	int overlapaa;
3115	double score, sumscore;
3116	int qstart, qend, tstart, tend;
3117	static char qal[N], tal[N], al[N];
3118	int nlocalhom;
3119
3120
3121
3122	count = 0;
3123	sumscore = 0.0;
3124	score = 0.0;
3125	nlocalhom = 0;
3126	while( 1 )
3127	{
3128
3129	if( feof( fp ) ) break;
3130
3131	while( fgets( b, B-1, fp ) )
3132	{
3133	if( !strncmp( " <Hit_def>", b, 19 ) \|\| !strncmp( " <Hsp_num>", b, 23 ) ) break;
3134	}
3135
3136	if( !strncmp( " <Hit_def>", b, 19 ) )
3137	{
3138	junban[count] = atoi( b+31 );
3139	nlocalhom = 0;
3140	}
3141
3142
3143	while( fgets( b, B-1, fp ) )
3144	if( !strncmp( " <Hsp_score>", b, 25 ) ) break;
3145	pt = b + 25;
3146	score = atof( pt );
3147	sumscore += score;
3148
3149
3150	while( fgets( b, B-1, fp ) )
3151	if( !strncmp( " <Hsp_query-from>", b, 30 ) ) break;
3152	pt = b + 30;
3153	qstart = atoi( pt ) - 1;
3154
3155
3156	while( fgets( b, B-1, fp ) )
3157	if( !strncmp( " <Hsp_query-to>", b, 28 ) ) break;
3158	pt = b + 28;
3159	qend = atoi( pt ) - 1;
3160
3161
3162	while( fgets( b, B-1, fp ) )
3163	if( !strncmp( " <Hsp_hit-from>", b, 28 ) ) break;
3164	pt = b + 28;
3165	tstart = atoi( pt ) - 1;
3166
3167
3168	while( fgets( b, B-1, fp ) )
3169	if( !strncmp( " <Hsp_hit-to>", b, 26 ) ) break;
3170	pt = b + 26;
3171	tend = atoi( pt ) - 1;
3172
3173
3174	while( fgets( b, B-1, fp ) )
3175	if( !strncmp( " <Hsp_align-len>", b, 29 ) ) break;
3176	pt = b + 29;
3177	overlapaa = atoi( pt );
3178
3179
3180	while( fgets( al, N-100, fp ) )
3181	if( !strncmp( " <Hsp_qseq>", al, 24 ) ) break;
3182
3183	strcpy( qal, al+24 );
3184	pt = qal;
3185	while( *++pt != '<' )
3186	;
3187	*pt = 0;
3188
3189
3190	while( fgets( al, N-100, fp ) )
3191	if( !strncmp( " <Hsp_hseq>", al, 24 ) ) break;
3192
3193	strcpy( tal, al+24 );
3194	pt = tal;
3195	while( *++pt != '<' )
3196	;
3197	*pt = 0;
3198
3199
3200	// fprintf( stderr, "t=%d, score = %f, qstart=%d, qend=%d, tstart=%d, tend=%d, overlapaa=%d\n", junban[count], score, qstart, qend, tstart, tend, overlapaa );
3201
3202	nlocalhom += addlocalhom_r( qal, tal, localhomlist+junban[count], qstart, tstart, score, overlapaa, nlocalhom );
3203
3204	while( fgets( b, B-1, fp ) )
3205	if( !strncmp( " </Hsp>:", b, 18 ) ) break;
3206
3207
3208	fgets( b, B-1, fp );
3209
3210
3211	if( !strncmp( " </Hit_hsps>", b, 21 ) )
3212	{
3213	dis[junban[count++]] = sumscore;
3214	sumscore = 0.0;
3215	fgets( b, B-1, fp );
3216	fgets( b, B-1, fp );
3217	if( !strncmp( " </Iteration_hits>", b, 23 ) ) break;
3218	}
3219	}
3220	return count;
3221	}
3222
3223	int ReadFasta34noalign( FILE fp, double dis, int qmem, char *name, LocalHom localhomlist )
3224	{
3225	int count=0;
3226	char b[B];
3227	char *pt;
3228	static int junban[M];
3229	int opt;
3230	double z, bits;
3231
3232
3233	count = 0;
3234	#if 0
3235	for( i=0; i<10000000 && count<nseq; i++ )
3236	#else
3237	while( !feof( fp ) )
3238	#endif
3239	{
3240	fgets( b, B-1, fp );
3241	if( !strncmp( "+==========+", b, 12 ) )
3242	{
3243	junban[count] = atoi( b+12 );
3244
3245	pt = strchr( b, ')' ) + 1;
3246	sscanf( pt, "%d %lf %lf", &opt, &bits, &z );
3247	dis[junban[count]] = (double)opt;
3248	count++;
3249
3250	}
3251	}
3252
3253	return count;
3254	}
3255	int ReadFasta34m10_nuc( FILE fp, double dis, int qmem, char *name, LocalHom localhomlist )
3256	{
3257	int count=0;
3258	char b[B];
3259	char *pt;
3260	static int junban[M];
3261	int overlapaa;
3262	int opt, qstart, qend, tstart, tend;
3263	double z, bits;
3264	int qal_display_start, tal_display_start;
3265	static char qal[N], tal[N];
3266	char qal2, tal2;
3267	int c;
3268
3269
3270	count = 0;
3271	#if 0
3272	for( i=0; i<10000000 && count<nseq; i++ )
3273	#else
3274	while( !feof( fp ) )
3275	#endif
3276	{
3277	fgets( b, B-1, fp );
3278	if( !strncmp( "+==========+", b, 12 ) )
3279	{
3280	junban[count] = atoi( b+12 );
3281
3282	if( strchr( b, 'r' ) ) continue;
3283
3284	pt = strchr( b, ']' ) + 1;
3285	sscanf( pt, "%d %lf %lf", &opt, &bits, &z );
3286	dis[junban[count]] = (double)opt;
3287	count++;
3288
3289	}
3290	else if( 0 == strncmp( ">>+==========+", b, 14 ) )
3291	{
3292	break;
3293	}
3294
3295	}
3296	if( !count ) return -1;
3297
3298	count = 0;
3299	while( 1 )
3300	{
3301	if( strncmp( ">>+==========+", b, 14 ) )
3302	{
3303	fgets( b, B-1, fp );
3304	if( feof( fp ) ) break;
3305	continue;
3306	}
3307	junban[count++] = atoi( b+14 );
3308	// fprintf( stderr, "t = %d\n", atoi( b+14 ) );
3309	while( fgets( b, B-1, fp ) )
3310	if( !strncmp( "; fa_opt:", b, 9 ) \|\| !strncmp( "; sw_s-w opt:", b, 13 ) ) break;
3311	pt = strstr( b, ":" ) +1;
3312	opt = atoi( pt );
3313
3314
3315	while( fgets( b, B-1, fp ) )
3316	if( !strncmp( "_overlap:", b+4, 9 ) ) break;
3317	pt = strstr( b, ":" ) +1;
3318	overlapaa = atoi( pt );
3319
3320	while( fgets( b, B-1, fp ) )
3321	if( !strncmp( "_start:", b+4, 7 ) ) break;
3322	pt = strstr( b, ":" ) +1;
3323	qstart = atoi( pt ) - 1;
3324
3325	while( fgets( b, B-1, fp ) )
3326	if( !strncmp( "_stop:", b+4, 6 ) ) break;
3327	pt = strstr( b, ":" ) +1;
3328	qend = atoi( pt ) - 1;
3329
3330	while( fgets( b, B-1, fp ) )
3331	if( !strncmp( "_display_start:", b+4, 15 ) ) break;
3332	pt = strstr( b, ":" ) +1;
3333	qal_display_start = atoi( pt ) - 1;
3334
3335	pt = qal;
3336	while( (c = fgetc( fp )) )
3337	{
3338	if( c == '>' )
3339	{
3340	ungetc( c, fp );
3341	break;
3342	}
3343	if( isalpha( c ) \|\| c == '-' )
3344	*pt++ = c;
3345	}
3346	*pt = 0;
3347
3348	while( fgets( b, B-1, fp ) )
3349	if( !strncmp( "_start:", b+4, 7 ) ) break;
3350	pt = strstr( b, ":" ) + 1;
3351	tstart = atoi( pt ) - 1;
3352
3353	while( fgets( b, B-1, fp ) )
3354	if( !strncmp( "_stop:", b+4, 6 ) ) break;
3355	pt = strstr( b, ":" ) + 1;
3356	tend = atoi( pt ) - 1;
3357
3358	while( fgets( b, B-1, fp ) )
3359	if( !strncmp( "_display_start:", b+4, 15 ) ) break;
3360	pt = strstr( b, ":" ) + 1;
3361	tal_display_start = atoi( pt ) - 1;
3362
3363	pt = tal;
3364	while( ( c = fgetc( fp ) ) )
3365	{
3366	if( c == '>' )
3367	{
3368	ungetc( c, fp );
3369	break;
3370	}
3371	if( isalpha( c ) \|\| c == '-' )
3372	*pt++ = c;
3373	}
3374	*pt = 0;
3375
3376	// fprintf( stderr, "(%d-%d:%d-%d)\n", qstart, qend, tstart, tend );
3377	// fprintf( stderr, "qal_display_start = %d, tal_display_start = %d\n", qal_display_start, tal_display_start );
3378
3379	// fprintf( stderr, "qal = %s\n", qal );
3380	// fprintf( stderr, "tal = %s\n", tal );
3381
3382	qal2 = cutal( qal, qal_display_start, qstart, qend );
3383	tal2 = cutal( tal, tal_display_start, tstart, tend );
3384
3385	// fprintf( stderr, "qal2 = %s\n", qal2 );
3386	// fprintf( stderr, "tal2 = %s\n", tal2 );
3387
3388	// fprintf( stderr, "putting %d - %d, opt = %d\n", qmem, junban[count-1], opt );
3389	putlocalhom( qal2, tal2, localhomlist+junban[count-1], qstart, tstart, opt, overlapaa );
3390	}
3391	// fprintf( stderr, "count = %d\n", count );
3392	return count;
3393	}
3394	int ReadFasta34m10( FILE fp, double dis, int qmem, char *name, LocalHom localhomlist )
3395	{
3396	int count=0;
3397	char b[B];
3398	char *pt;
3399	static int junban[M];
3400	int overlapaa;
3401	int opt, qstart, qend, tstart, tend;
3402	double z, bits;
3403	int qal_display_start, tal_display_start;
3404	static char qal[N], tal[N];
3405	char qal2, tal2;
3406	int c;
3407
3408
3409	count = 0;
3410	#if 0
3411	for( i=0; i<10000000 && count<nseq; i++ )
3412	#else
3413	while( !feof( fp ) )
3414	#endif
3415	{
3416	fgets( b, B-1, fp );
3417	if( !strncmp( "+==========+", b, 12 ) )
3418	{
3419	junban[count] = atoi( b+12 );
3420
3421	pt = strchr( b, ')' ) + 1;
3422	sscanf( pt, "%d %lf %lf", &opt, &bits, &z );
3423	dis[junban[count]] = (double)opt;
3424	count++;
3425
3426	}
3427	else if( 0 == strncmp( ">>+==========+", b, 14 ) )
3428	{
3429	break;
3430	}
3431
3432	}
3433	if( !count ) return -1;
3434
3435	count = 0;
3436	while( 1 )
3437	{
3438	if( strncmp( ">>+==========+", b, 14 ) )
3439	{
3440	fgets( b, B-1, fp );
3441	if( feof( fp ) ) break;
3442	continue;
3443	}
3444	junban[count++] = atoi( b+14 );
3445	// fprintf( stderr, "t = %d\n", atoi( b+14 ) );
3446	while( fgets( b, B-1, fp ) )
3447	if( !strncmp( "; fa_opt:", b, 9 ) \|\| !strncmp( "; sw_s-w opt:", b, 13 ) ) break;
3448	pt = strstr( b, ":" ) +1;
3449	opt = atoi( pt );
3450
3451
3452	while( fgets( b, B-1, fp ) )
3453	if( !strncmp( "_overlap:", b+4, 9 ) ) break;
3454	pt = strstr( b, ":" ) +1;
3455	overlapaa = atoi( pt );
3456
3457	while( fgets( b, B-1, fp ) )
3458	if( !strncmp( "_start:", b+4, 7 ) ) break;
3459	pt = strstr( b, ":" ) +1;
3460	qstart = atoi( pt ) - 1;
3461
3462	while( fgets( b, B-1, fp ) )
3463	if( !strncmp( "_stop:", b+4, 6 ) ) break;
3464	pt = strstr( b, ":" ) +1;
3465	qend = atoi( pt ) - 1;
3466
3467	while( fgets( b, B-1, fp ) )
3468	if( !strncmp( "_display_start:", b+4, 15 ) ) break;
3469	pt = strstr( b, ":" ) +1;
3470	qal_display_start = atoi( pt ) - 1;
3471
3472	pt = qal;
3473	while( (c = fgetc( fp )) )
3474	{
3475	if( c == '>' )
3476	{
3477	ungetc( c, fp );
3478	break;
3479	}
3480	if( isalpha( c ) \|\| c == '-' )
3481	*pt++ = c;
3482	}
3483	*pt = 0;
3484
3485	while( fgets( b, B-1, fp ) )
3486	if( !strncmp( "_start:", b+4, 7 ) ) break;
3487	pt = strstr( b, ":" ) + 1;
3488	tstart = atoi( pt ) - 1;
3489
3490	while( fgets( b, B-1, fp ) )
3491	if( !strncmp( "_stop:", b+4, 6 ) ) break;
3492	pt = strstr( b, ":" ) + 1;
3493	tend = atoi( pt ) - 1;
3494
3495	while( fgets( b, B-1, fp ) )
3496	if( !strncmp( "_display_start:", b+4, 15 ) ) break;
3497	pt = strstr( b, ":" ) + 1;
3498	tal_display_start = atoi( pt ) - 1;
3499
3500	pt = tal;
3501	while( ( c = fgetc( fp ) ) )
3502	{
3503	if( c == '>' )
3504	{
3505	ungetc( c, fp );
3506	break;
3507	}
3508	if( isalpha( c ) \|\| c == '-' )
3509	*pt++ = c;
3510	}
3511	*pt = 0;
3512
3513	// fprintf( stderr, "(%d-%d:%d-%d)\n", qstart, qend, tstart, tend );
3514	// fprintf( stderr, "qal_display_start = %d, tal_display_start = %d\n", qal_display_start, tal_display_start );
3515
3516	// fprintf( stderr, "qal = %s\n", qal );
3517	// fprintf( stderr, "tal = %s\n", tal );
3518
3519	qal2 = cutal( qal, qal_display_start, qstart, qend );
3520	tal2 = cutal( tal, tal_display_start, tstart, tend );
3521
3522	// fprintf( stderr, "qal2 = %s\n", qal2 );
3523	// fprintf( stderr, "tal2 = %s\n", tal2 );
3524
3525	// fprintf( stderr, "putting %d - %d, opt = %d\n", qmem, junban[count-1], opt );
3526	putlocalhom( qal2, tal2, localhomlist+junban[count-1], qstart, tstart, opt, overlapaa );
3527	}
3528	// fprintf( stderr, "count = %d\n", count );
3529	return count;
3530	}
3531	int ReadFasta34m10_scoreonly_nucbk( FILE fp, double dis, int nin )
3532	{
3533	int count=0;
3534	char b[B];
3535	char *pt;
3536	int pos;
3537	int opt;
3538	double z, bits;
3539
3540	count = 0;
3541	while( !feof( fp ) )
3542	{
3543	fgets( b, B-1, fp );
3544	if( !strncmp( "+===========+", b, 13 ) )
3545	{
3546	pos = atoi( b+13 );
3547
3548	if( strchr( b, 'r' ) ) continue;
3549
3550	// pt = strchr( b, ')' ) + 1;
3551	pt = strchr( b, ']' ) + 1;
3552	sscanf( pt, "%d %lf %lf", &opt, &bits, &z );
3553	dis[pos] += (double)opt;
3554	count++;
3555	#if 0
3556	fprintf( stderr, "b=%s\n", b );
3557	fprintf( stderr, "opt=%d\n", opt );
3558	fprintf( stderr, "pos=%d\n", pos );
3559	fprintf( stderr, "dis[pos]=%f\n", dis[pos] );
3560	#endif
3561
3562	}
3563	else if( 0 == strncmp( ">>><<<", b, 6 ) )
3564	{
3565	break;
3566	}
3567
3568	}
3569	if( !count ) return -1;
3570
3571	return count;
3572	}
3573
3574	int ReadFasta34m10_scoreonly_nuc( FILE fp, double dis, int nin )
3575	{
3576	int count=0;
3577	char b[B];
3578	char *pt;
3579	int pos;
3580	int opt;
3581	double z, bits;
3582	int c;
3583	int *yonda;
3584
3585
3586	yonda = AllocateIntVec( nin );
3587	for( c=0; c<nin; c++ ) yonda[c] = 0;
3588	for( c=0; c<nin; c++ ) dis[c] = 0.0;
3589
3590	count = 0;
3591	while( !feof( fp ) )
3592	{
3593	fgets( b, B-1, fp );
3594	if( !strncmp( "+===========+", b, 13 ) )
3595	{
3596	pos = atoi( b+13 );
3597
3598	if( strchr( b, 'r' ) ) continue;
3599
3600	// pt = strchr( b, ')' ) + 1;
3601	pt = strchr( b, ']' ) + 1;
3602	sscanf( pt, "%d %lf %lf", &opt, &bits, &z );
3603	if( yonda[pos] == 0 )
3604	{
3605	dis[pos] += (double)opt;
3606	yonda[pos] = 1;
3607	}
3608	count++;
3609	#if 0
3610	fprintf( stderr, "b=%s\n", b );
3611	fprintf( stderr, "opt=%d\n", opt );
3612	fprintf( stderr, "pos=%d\n", pos );
3613	fprintf( stderr, "dis[pos]=%f\n", dis[pos] );
3614	#endif
3615
3616	}
3617	else if( !strncmp( ">>>", b, 3 ) )
3618	{
3619	for( c=0; c<nin; c++ ) yonda[c] = 0;
3620	}
3621	else if( 0 == strncmp( ">>><<<", b, 6 ) )
3622	{
3623	break;
3624	}
3625
3626	}
3627
3628	free( yonda );
3629
3630	if( !count ) return -1;
3631
3632	return count;
3633	}
3634
3635	int ReadFasta34m10_scoreonly( FILE fp, double dis, int nin )
3636	{
3637	int count=0;
3638	char b[B];
3639	char *pt;
3640	int pos;
3641	int opt;
3642	double z, bits;
3643	int c;
3644	int *yonda;
3645
3646
3647	yonda = AllocateIntVec( nin );
3648	for( c=0; c<nin; c++ ) yonda[c] = 0;
3649	for( c=0; c<nin; c++ ) dis[c] = 0.0;
3650
3651	count = 0;
3652	while( !feof( fp ) )
3653	{
3654	fgets( b, B-1, fp );
3655	if( !strncmp( "+===========+", b, 13 ) )
3656	{
3657	pos = atoi( b+13 );
3658
3659	pt = strchr( b, ')' ) + 1;
3660	sscanf( pt, "%d %lf %lf", &opt, &bits, &z );
3661	if( yonda[pos] == 0 )
3662	{
3663	dis[pos] += (double)opt;
3664	yonda[pos] = 1;
3665	}
3666	count++;
3667	#if 0
3668	fprintf( stderr, "b=%s\n", b );
3669	fprintf( stderr, "opt=%d\n", opt );
3670	fprintf( stderr, "pos=%d\n", pos );
3671	fprintf( stderr, "dis[pos]=%f\n", dis[pos] );
3672	#endif
3673
3674	}
3675	else if( !strncmp( ">>>", b, 3 ) )
3676	{
3677	for( c=0; c<nin; c++ ) yonda[c] = 0;
3678	}
3679	else if( 0 == strncmp( ">>><<<", b, 6 ) )
3680	{
3681	break;
3682	}
3683
3684	}
3685
3686	free( yonda );
3687
3688	if( !count ) return -1;
3689
3690	return count;
3691	}
3692	int ReadFasta34( FILE fp, double dis, int nseq, char name[M][B], LocalHom *localhomlist )
3693	{
3694	int count=0;
3695	char b[B];
3696	char *pt;
3697	static int junban[M];
3698	int overlapaa;
3699	int opt, qstart, qend, tstart, tend;
3700	double z, bits;
3701
3702
3703	count = 0;
3704	#if 0
3705	for( i=0; i<10000000 && count<nseq; i++ )
3706	#else
3707	while( !feof( fp ) )
3708	#endif
3709	{
3710	fgets( b, B-1, fp );
3711	if( !strncmp( "+==========+", b, 12 ) )
3712	{
3713	junban[count] = atoi( b+12 );
3714
3715	pt = strchr( b, ')' ) + 1;
3716	sscanf( pt, "%d %lf %lf", &opt, &bits, &z );
3717	dis[junban[count]] = (double)opt;
3718	count++;
3719
3720	}
3721	else if( 0 == strncmp( ">>+==========+", b, 14 ) )
3722	{
3723	break;
3724	}
3725
3726	}
3727	if( !count ) return -1;
3728
3729	count = 0;
3730	while( !feof( fp ) )
3731	{
3732	if( !strncmp(">>+==========+", b, 14 ) )
3733	{
3734	junban[count] = atoi( b+14 );
3735	count++;
3736	fgets( b, B-1, fp ); // initn:
3737	pt = strstr( b, "opt: " ) + 5;
3738	localhomlist[junban[count-1]].opt = atof( pt );
3739	fgets( b, B-1, fp ); // Smith-Waterman score
3740	pt = strstr( b, "ungapped) in " ) + 13;
3741	sscanf( pt, "%d", &overlapaa );
3742	fprintf( stderr, "pt = %s, overlapaa = %d\n", pt, overlapaa );
3743	pt = strstr( b, "overlap (" ) + 8;
3744	sscanf( pt, "(%d-%d:%d-%d)", &qstart, &qend, &tstart, &tend );
3745	localhomlist[junban[count-1]].overlapaa = overlapaa;
3746	localhomlist[junban[count-1]].start1 = qstart-1;
3747	localhomlist[junban[count-1]].end1 = qend-1;
3748	localhomlist[junban[count-1]].start2 = tstart-1;
3749	localhomlist[junban[count-1]].end2 = tend-1;
3750	}
3751	fgets( b, B-1, fp );
3752	}
3753	fprintf( stderr, "count = %d\n", count );
3754	return count;
3755	}
3756
3757	int ReadFasta3( FILE fp, double dis, int nseq, char name[M][B] )
3758	{
3759	int count=0;
3760	char b[B];
3761	char *pt;
3762	int junban[M];
3763	int initn, init1, opt;
3764	double z;
3765
3766	count = 0;
3767	#if 0
3768	for( i=0; i<10000000 && count<nseq; i++ )
3769	#else
3770	while( !feof( fp ) )
3771	#endif
3772	{
3773	fgets( b, B-1, fp );
3774	if( !strncmp( "+==========+", b, 12 ) )
3775	{
3776	junban[count] = atoi( b+12 );
3777
3778	pt = strchr( b, ')' ) + 1;
3779	sscanf( pt, "%d %d %d %lf", &initn, &init1, &opt, &z );
3780	dis[junban[count]] = (double)opt;
3781	count++;
3782	}
3783	}
3784	return 0;
3785	}
3786
3787	int ReadFasta( FILE fp, double dis, int nseq, char name[M][B] )
3788	{
3789	int i, count=0;
3790	char b[B];
3791	int junban[M];
3792	int initn, init1, opt;
3793
3794	count = 0;
3795	for( i=0; i<nseq; i++ ) dis[i] = 0.0;
3796	for( i=0; !feof( fp ) && count<nseq; i++ )
3797	{
3798	fgets( b, B-1, fp );
3799	if( !strncmp( "+==========+", b, 12 ) )
3800	{
3801	junban[count] = atoi( b+12 );
3802
3803	sscanf( b+50, "%d %d %d", &initn, &init1, &opt );
3804	dis[junban[count]] = (double)opt;
3805	count++;
3806	}
3807	}
3808
3809	/*
3810	count = 0;
3811	for( i=0; i<100000 && count<nseq; i++ )
3812	{
3813	fgets( b, B-1, fp );
3814	if( !strncmp( name[junban[count]], b, 20 ) )
3815	{
3816	dis[junban[count]] = atof( b+65 );
3817	count++;
3818	}
3819	}
3820	*/
3821	return 0;
3822	}
3823
3824
3825	int ReadOpt( FILE *fp, int opt[M], int nseq, char name[M][B] )
3826	{
3827	int i, count=0;
3828	char b[B];
3829	int junban[M];
3830	int optt, initn, init1;
3831
3832	count = 0;
3833	for( i=0; i<10000000 && count<nseq; i++ )
3834	{
3835	fgets( b, B-1, fp );
3836	if( !strncmp( "+==========+", b, 12 ) )
3837	{
3838	junban[count] = atoi( b+12 );
3839	sscanf( b+50, "%d %d %d", &initn, &init1, &optt );
3840	opt[junban[count]] = (double)optt;
3841	count++;
3842	}
3843	}
3844	return 0;
3845	}
3846
3847	int ReadOpt2( FILE *fp, int opt[M], int nseq, char name[M][B] )
3848	{
3849	int i, count=0;
3850	char b[B];
3851	int junban[M];
3852
3853	count = 0;
3854	for( i=0; i<10000000 && count<nseq; i++ )
3855	{
3856	fgets( b, B-1, fp );
3857	if( !strncmp( "+==========+", b, 12 ) )
3858	{
3859	junban[count] = atoi( b+12 );
3860	opt[junban[count]] = atoi( b+65 );
3861	count++;
3862	}
3863	}
3864	return 0;
3865	}
3866
3867
3868
3869	int writePre( int nseq, char name, int nlen[M], char aseq, int force )
3870	{
3871	#if USE_XCED
3872	int i, value;
3873	if( !signalSM )
3874	{
3875	if( force )
3876	{
3877	rewind( prep_g );
3878	if( devide ) dvWrite( prep_g, nseq, name, nlen, aseq );
3879	#if 0
3880	else WriteGapFill( prep_g, nseq, name, nlen, aseq );
3881	#else
3882	else writeData( prep_g, nseq, name, nlen, aseq );
3883	#endif
3884	}
3885	return( 0 );
3886	}
3887	for( i=0; i<10; i++ )
3888	{
3889	#if IODEBUG
3890	fprintf( stderr, "SEMAPHORE = %d\n", signalSM[SEMAPHORE] );
3891	#endif
3892	if( signalSM[SEMAPHORE]-- > 0 )
3893	{
3894	#if 0 /* /tmp/pre €ÎŽØ·ž€Ç€Ï€º€·€¿ */
3895	if( ferror( prep_g ) ) prep_g = fopen( "pre", "w" );
3896	if( !prep_g ) ErrorExit( "Cannot re-open pre." );
3897	#endif
3898	rewind( prep_g );
3899	signalSM[STATUS] = IMA_KAITERU;
3900	#if IODEBUG
3901	if( force ) fprintf( stderr, "FINAL " );
3902	#endif
3903	if( devide ) dvWrite( prep_g, nseq, name, nlen, aseq );
3904	else WriteGapFill( prep_g, nseq, name, nlen, aseq );
3905	/*
3906	fprintf( prep_g, '\EOF' );
3907	*/
3908	fflush( prep_g );
3909	if( force ) signalSM[STATUS] = OSHIMAI;
3910	else signalSM[STATUS] = KAKIOWATTA;
3911	value = 1;
3912	signalSM[SEMAPHORE]++;
3913	#if IODEBUG
3914	fprintf( stderr, "signalSM[STATUS] = %c\n", signalSM[STATUS] );
3915	#endif
3916	break;
3917	}
3918	else
3919	{
3920	#if IODEBUG
3921	fprintf( stderr, "YONDERUKARA_AKIRAMERU\n" );
3922	#endif
3923	value = 0;
3924	signalSM[SEMAPHORE]++;
3925	if( !force ) break;
3926	#if IODEBUG
3927	fprintf( stderr, "MATSU\n" );
3928	#endif
3929	sleep( 1 );
3930	}
3931	}
3932	if( force && !value ) ErrorExit( "xced ga pre wo hanasanai \n" );
3933	return( value );
3934	#else
3935	if( force )
3936	{
3937	rewind( prep_g );
3938	writeData_pointer( prep_g, nseq, name, nlen, aseq );
3939	}
3940	#endif
3941	return( 0 );
3942	}
3943
3944
3945	void readOtherOptions( int ppidptr, int fftThresholdptr, int *fftWinSizeptr )
3946	{
3947	if( calledByXced )
3948	{
3949	FILE *fp = fopen( "pre", "r" );
3950	char b[B];
3951	if( !fp ) ErrorExit( "Cannot open pre.\n" );
3952	fgets( b, B-1, fp );
3953	sscanf( b, "%d %d %d", ppidptr, fftThresholdptr, fftWinSizeptr );
3954	fclose( fp );
3955	#if IODEBUG
3956	fprintf( stderr, "b = %s\n", b );
3957	fprintf( stderr, "ppid = %d\n", ppid );
3958	fprintf( stderr, "fftThreshold = %d\n", fftThreshold );
3959	fprintf( stderr, "fftWinSize = %d\n", fftWinSize );
3960	#endif
3961	}
3962	else
3963	{
3964	*ppidptr = 0;
3965	*fftThresholdptr = FFT_THRESHOLD;
3966	if( dorp == 'd' )
3967	*fftWinSizeptr = FFT_WINSIZE_D;
3968	else
3969	*fftWinSizeptr = FFT_WINSIZE_P;
3970	}
3971	#if 0
3972	fprintf( stderr, "fftThresholdptr=%d\n", *fftThresholdptr );
3973	fprintf( stderr, "fftWinSizeptr=%d\n", *fftWinSizeptr );
3974	#endif
3975	}
3976
3977	void initSignalSM( void )
3978	{
3979	// int signalsmid;
3980
3981	#if IODEBUG
3982	if( ppid ) fprintf( stderr, "PID of xced = %d\n", ppid );
3983	#endif
3984	if( !ppid )
3985	{
3986	signalSM = NULL;
3987	return;
3988	}
3989
3990	#if 0
3991	signalsmid = shmget( (key_t)ppid, 3, IPC_ALLOC \| 0666 );
3992	if( signalsmid == -1 ) ErrorExit( "Cannot get Shared memory for signal.\n" );
3993	signalSM = shmat( signalsmid, 0, 0 );
3994	if( (int)signalSM == -1 ) ErrorExit( "Cannot attatch Shared Memory for signal!\n" );
3995	signalSM[STATUS] = IMA_KAITERU;
3996	signalSM[SEMAPHORE] = 1;
3997	#endif
3998	}
3999
4000	void initFiles( void )
4001	{
4002	char pname[100];
4003	if( ppid )
4004	sprintf( pname, "/tmp/pre.%d", ppid );
4005	else
4006	sprintf( pname, "pre" );
4007	prep_g = fopen( pname, "w" );
4008	if( !prep_g ) ErrorExit( "Cannot open pre" );
4009
4010	trap_g = fopen( "trace", "w" );
4011	if( !trap_g ) ErrorExit( "cannot open trace" );
4012	fprintf( trap_g, "PID = %d\n", getpid() );
4013	fflush( trap_g );
4014	}
4015
4016
4017	void WriteForFasta( FILE fp, int locnjob, char name, int nlen[M], char *aseq )
4018	{
4019	static char b[N];
4020	int i, j;
4021	int nalen[M];
4022
4023	for( i=0; i<locnjob; i++ )
4024	{
4025	nalen[i] = strlen( aseq[i] );
4026	fprintf( fp, ">%s\n", name[i] );
4027	for( j=0; j<nalen[i]; j=j+C )
4028	{
4029	strncpy( b, aseq[i]+j, C ); b[C] = 0;
4030	fprintf( fp, "%s\n",b );
4031	}
4032	}
4033	}
4034
4035	void readlocalhomtable2( FILEfp, int njob, LocalHom localhomtable, char kozoarivec )
4036	{
4037	double opt;
4038	static char buff[B];
4039	char infor[100];
4040	int i, j, overlapaa, start1, end1, start2, end2;
4041	LocalHom tmpptr1, tmpptr2;
4042
4043	// for( i=0; i<njob; i++ ) for( j=0; j<njob; j++ ) nlocalhom[i][j] = 0;
4044
4045	while ( NULL != fgets( buff, B-1, fp ) )
4046	{
4047	// fprintf( stderr, "\n" );
4048	sscanf( buff, "%d %d %d %lf %d %d %d %d %s", &i, &j, &overlapaa, &opt, &start1, &end1, &start2, &end2, infor );
4049	if( *infor == 'k' ) kozoarivec[i] = kozoarivec[j] = 1;
4050
4051	#if 0
4052	if( start1 == end1 \|\| start2 == end2 ) continue; //mondai ari
4053	#endif
4054
4055	// if( i < j )
4056	{
4057	if( localhomtable[i][j].nokori++ > 0 )
4058	{
4059	tmpptr1 = localhomtable[i][j].last;
4060	// fprintf( stderr, "reallocating, localhomtable[%d][%d].nokori = %d\n", i, j, localhomtable[i][j].nokori );
4061	tmpptr1->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
4062	tmpptr1 = tmpptr1->next;
4063	tmpptr1->extended = -1;
4064	tmpptr1->next = NULL;
4065	localhomtable[i][j].last = tmpptr1;
4066	// fprintf( stderr, "### i,j = %d,%d, nokori=%d\n", i, j, localhomtable[i][j].nokori );
4067	}
4068	else
4069	{
4070	tmpptr1 = localhomtable[i]+j;
4071	// fprintf( stderr, "### i,j = %d,%d, nokori=%d\n", i, j, localhomtable[i][j].nokori );
4072	}
4073
4074	tmpptr1->start1 = start1;
4075	tmpptr1->start2 = start2;
4076	tmpptr1->end1 = end1;
4077	tmpptr1->end2 = end2;
4078	// tmpptr1->opt = ( opt / overlapaa + 0.00 ) / 5.8 * 600;
4079	// tmpptr1->opt = opt;
4080	tmpptr1->opt = ( opt + 0.00 ) / 5.8 * 600;
4081	tmpptr1->overlapaa = overlapaa;
4082	tmpptr1->korh = *infor;
4083	}
4084	// else
4085	{
4086	if( localhomtable[j][i].nokori++ > 0 )
4087	{
4088	tmpptr2 = localhomtable[j][i].last;
4089	tmpptr2->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
4090	tmpptr2 = tmpptr2->next;
4091	tmpptr2->extended = -1;
4092	tmpptr2->next = NULL;
4093	localhomtable[j][i].last = tmpptr2;
4094	// fprintf( stderr, "### i,j = %d,%d, nokori=%d\n", j, i, localhomtable[j][i].nokori );
4095	}
4096	else
4097	{
4098	tmpptr2 = localhomtable[j]+i;
4099	// fprintf( stderr, "### i,j = %d,%d, nokori=%d\n", j, i, localhomtable[j][i].nokori );
4100	}
4101
4102	tmpptr2->start2 = start1;
4103	tmpptr2->start1 = start2;
4104	tmpptr2->end2 = end1;
4105	tmpptr2->end1 = end2;
4106	// tmpptr2->opt = ( opt / overlapaa + 0.00 ) / 5.8 * 600;
4107	// tmpptr2->opt = opt;
4108	tmpptr2->opt = ( opt + 0.00 ) / 5.8 * 600;
4109	tmpptr2->overlapaa = overlapaa;
4110	tmpptr2->korh = *infor;
4111
4112	// fprintf( stderr, "i=%d, j=%d, st1=%d, en1=%d, opt = %f\n", i, j, tmpptr1->start1, tmpptr1->end1, opt );
4113	}
4114
4115	}
4116	}
4117	void readlocalhomtable( FILEfp, int njob, LocalHom localhomtable, char kozoarivec )
4118	{
4119	double opt;
4120	static char buff[B];
4121	char infor[100];
4122	int i, j, overlapaa, start1, end1, start2, end2;
4123	int **nlocalhom = NULL;
4124	LocalHom tmpptr1=NULL, tmpptr2=NULL; // by D.Mathog, a guess
4125
4126	nlocalhom = AllocateIntMtx( njob, njob );
4127	for( i=0; i<njob; i++ ) for( j=0; j<njob; j++ ) nlocalhom[i][j] = 0;
4128
4129	while ( NULL != fgets( buff, B-1, fp ) )
4130	{
4131	// fprintf( stderr, "\n" );
4132	sscanf( buff, "%d %d %d %lf %d %d %d %d %s", &i, &j, &overlapaa, &opt, &start1, &end1, &start2, &end2, infor );
4133	if( *infor == 'k' ) kozoarivec[i] = kozoarivec[j] = 1;
4134
4135	#if 0
4136	if( start1 == end1 \|\| start2 == end2 ) continue; //mondai ari
4137	#endif
4138
4139
4140	// if( i < j )
4141	{
4142	if( nlocalhom[i][j]++ > 0 )
4143	{
4144	// fprintf( stderr, "reallocating, nlocalhom[%d][%d] = %d\n", i, j, nlocalhom[i][j] );
4145	tmpptr1->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
4146	tmpptr1 = tmpptr1->next;
4147	tmpptr1->next = NULL;
4148	}
4149	else
4150	{
4151	tmpptr1 = localhomtable[i]+j;
4152	// fprintf( stderr, "nlocalhom[%d][%d] = %d\n", i, j, nlocalhom[i][j] );
4153	}
4154
4155	tmpptr1->start1 = start1; // CHUUI!!!!
4156	tmpptr1->start2 = start2;
4157	tmpptr1->end1 = end1; // CHUUI!!!!
4158	tmpptr1->end2 = end2;
4159	// tmpptr1->opt = ( opt / overlapaa + 0.00 ) / 5.8 * 600;
4160	// tmpptr1->opt = opt;
4161	tmpptr1->opt = ( opt + 0.00 ) / 5.8 * 600;
4162	tmpptr1->overlapaa = overlapaa;
4163	tmpptr1->korh = *infor;
4164
4165	// fprintf( stderr, "i=%d, j=%d, opt = %f\n", i, j, opt );
4166	}
4167	// else
4168	{
4169	if( nlocalhom[j][i]++ > 0 )
4170	{
4171	tmpptr2->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
4172	tmpptr2 = tmpptr2->next;
4173	tmpptr2->next = NULL;
4174	}
4175	else
4176	tmpptr2 = localhomtable[j]+i;
4177
4178	tmpptr2->start2 = start1; // CHUUI!!!!
4179	tmpptr2->start1 = start2;
4180	tmpptr2->end2 = end1; // CHUUI!!!!
4181	tmpptr2->end1 = end2;
4182	// tmpptr2->opt = ( opt / overlapaa + 0.00 ) / 5.8 * 600;
4183	// tmpptr2->opt = opt;
4184	tmpptr2->opt = ( opt + 0.00 ) / 5.8 * 600;
4185	tmpptr2->overlapaa = overlapaa;
4186	tmpptr2->korh = *infor;
4187
4188	// fprintf( stderr, "j=%d, i=%d, opt = %f\n", j, i, opt );
4189	}
4190
4191	}
4192	FreeIntMtx( nlocalhom );
4193	}
4194
4195
4196	void readlocalhomtable_two( FILEfp, int norg, int nadd, LocalHom localhomtable, LocalHom localhomtablex, char kozoarivec ) // for test only
4197	{
4198	double opt;
4199	static char buff[B];
4200	char infor[100];
4201	int i, j, overlapaa, start1, end1, start2, end2;
4202	int **nlocalhom = NULL;
4203	int **nlocalhomx = NULL;
4204	LocalHom tmpptr1=NULL, tmpptr2=NULL; // by D.Mathog, a guess
4205
4206	nlocalhom = AllocateIntMtx( norg, nadd );
4207	for( i=0; i<norg; i++ ) for( j=0; j<nadd; j++ ) nlocalhom[i][j] = 0;
4208	nlocalhomx = AllocateIntMtx( nadd, norg );
4209	for( i=0; i<nadd; i++ ) for( j=0; j<norg; j++ ) nlocalhomx[i][j] = 0;
4210
4211	while ( NULL != fgets( buff, B-1, fp ) )
4212	{
4213	// fprintf( stderr, "\n" );
4214	sscanf( buff, "%d %d %d %lf %d %d %d %d %s", &i, &j, &overlapaa, &opt, &start1, &end1, &start2, &end2, infor );
4215	if( *infor == 'k' )
4216	{
4217	fprintf( stderr, "Not supported!\n" );
4218	exit( 1 );
4219	}
4220	j -= norg;
4221
4222	#if 0
4223	if( start1 == end1 \|\| start2 == end2 ) continue; //mondai ari
4224	#endif
4225
4226
4227	// if( i < j )
4228	{
4229	if( nlocalhom[i][j]++ > 0 )
4230	{
4231	// fprintf( stderr, "reallocating, nlocalhom[%d][%d] = %d\n", i, j, nlocalhom[i][j] );
4232	tmpptr1->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
4233	tmpptr1 = tmpptr1->next;
4234	tmpptr1->next = NULL;
4235	}
4236	else
4237	{
4238	tmpptr1 = localhomtable[i]+j;
4239	// fprintf( stderr, "nlocalhom[%d][%d] = %d\n", i, j, nlocalhom[i][j] );
4240	}
4241
4242	tmpptr1->start1 = start1; // CHUUI!!!!
4243	tmpptr1->start2 = start2;
4244	tmpptr1->end1 = end1; // CHUUI!!!!
4245	tmpptr1->end2 = end2;
4246	// tmpptr1->opt = ( opt / overlapaa + 0.00 ) / 5.8 * 600;
4247	// tmpptr1->opt = opt;
4248	tmpptr1->opt = ( opt + 0.00 ) / 5.8 * 600;
4249	tmpptr1->overlapaa = overlapaa;
4250	tmpptr1->korh = *infor;
4251
4252	// fprintf( stderr, "i=%d, j=%d, opt = %f\n", i, j, opt );
4253	}
4254
4255	{
4256	if( nlocalhomx[j][i]++ > 0 )
4257	{
4258	tmpptr2->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
4259	tmpptr2 = tmpptr2->next;
4260	tmpptr2->next = NULL;
4261	}
4262	else
4263	tmpptr2 = localhomtablex[j]+i;
4264
4265	tmpptr2->start2 = start1+1; // CHUUI!!!!
4266	tmpptr2->start1 = start2;
4267	tmpptr2->end2 = end1+1; // CHUUI!!!!
4268	tmpptr2->end1 = end2;
4269	// tmpptr2->opt = ( opt / overlapaa + 0.00 ) / 5.8 * 600;
4270	// tmpptr2->opt = opt;
4271	tmpptr2->opt = ( opt + 0.00 ) / 5.8 * 600;
4272	tmpptr2->overlapaa = overlapaa;
4273	tmpptr2->korh = *infor;
4274
4275	// fprintf( stderr, "j=%d, i=%d, opt = %f\n", j, i, opt );
4276	}
4277
4278	}
4279	FreeIntMtx( nlocalhom );
4280	FreeIntMtx( nlocalhomx );
4281	}
4282
4283	void readlocalhomtable_one( FILEfp, int norg, int nadd, LocalHom localhomtable, char kozoarivec ) // for test only
4284	{
4285	double opt;
4286	static char buff[B];
4287	char infor[100];
4288	int i, j, overlapaa, start1, end1, start2, end2;
4289	int **nlocalhom = NULL;
4290	LocalHom *tmpptr1=NULL; // by D.Mathog, a guess
4291
4292	nlocalhom = AllocateIntMtx( norg, nadd );
4293	for( i=0; i<norg; i++ ) for( j=0; j<nadd; j++ ) nlocalhom[i][j] = 0;
4294
4295	while ( NULL != fgets( buff, B-1, fp ) )
4296	{
4297	// fprintf( stderr, "\n" );
4298	sscanf( buff, "%d %d %d %lf %d %d %d %d %s", &i, &j, &overlapaa, &opt, &start1, &end1, &start2, &end2, infor );
4299	if( *infor == 'k' )
4300	{
4301	fprintf( stderr, "Not supported!\n" );
4302	exit( 1 );
4303	}
4304	j -= norg;
4305
4306	#if 0
4307	if( start1 == end1 \|\| start2 == end2 ) continue; //mondai ari
4308	#endif
4309
4310
4311	// if( i < j )
4312	{
4313	if( nlocalhom[i][j]++ > 0 )
4314	{
4315	// fprintf( stderr, "reallocating, nlocalhom[%d][%d] = %d\n", i, j, nlocalhom[i][j] );
4316	tmpptr1->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
4317	tmpptr1 = tmpptr1->next;
4318	tmpptr1->next = NULL;
4319	}
4320	else
4321	{
4322	tmpptr1 = localhomtable[i]+j;
4323	// fprintf( stderr, "nlocalhom[%d][%d] = %d\n", i, j, nlocalhom[i][j] );
4324	}
4325
4326	tmpptr1->start1 = start1; // CHUUI!!!!
4327	tmpptr1->start2 = start2;
4328	tmpptr1->end1 = end1; // CHUUI!!!!
4329	tmpptr1->end2 = end2;
4330	// tmpptr1->opt = ( opt / overlapaa + 0.00 ) / 5.8 * 600;
4331	// tmpptr1->opt = opt;
4332	tmpptr1->opt = ( opt + 0.00 ) / 5.8 * 600;
4333	tmpptr1->overlapaa = overlapaa;
4334	tmpptr1->korh = *infor;
4335
4336	// fprintf( stderr, "i=%d, j=%d, opt = %f\n", i, j, opt );
4337	}
4338
4339	}
4340	FreeIntMtx( nlocalhom );
4341	}
4342
4343	void outlocalhom_part( LocalHom **localhom, int norg, int nadd )
4344	{
4345	int i, j;
4346	LocalHom *tmpptr;
4347	for( i=0; i<norg; i++ ) for( j=0; j<nadd; j++ )
4348	{
4349	tmpptr = localhom[i]+j;
4350	fprintf( stderr, "%d-%d\n", i, j+norg );
4351	do
4352	{
4353	fprintf( stderr, "reg1=%d-%d, reg2=%d-%d, imp=%f, opt=%f\n", tmpptr->start1, tmpptr->end1, tmpptr->start2, tmpptr->end2, tmpptr->importance, tmpptr->opt );
4354	}
4355	while( (tmpptr=tmpptr->next) );
4356	}
4357	}
4358
4359	void outlocalhom( LocalHom **localhom, int nseq )
4360	{
4361	int i, j;
4362	LocalHom *tmpptr;
4363	for( i=0; i<nseq; i++ ) for( j=0; j<nseq; j++ )
4364	{
4365	tmpptr = localhom[i]+j;
4366	fprintf( stderr, "%d-%d\n", i, j );
4367	do
4368	{
4369	fprintf( stderr, "reg1=%d-%d, reg2=%d-%d, imp=%f, opt=%f\n", tmpptr->start1, tmpptr->end1, tmpptr->start2, tmpptr->end2, tmpptr->importance, tmpptr->opt );
4370	}
4371	while( (tmpptr=tmpptr->next) );
4372	}
4373	}
4374
4375	void outlocalhompt( LocalHom ***localhom, int n1, int n2 )
4376	{
4377	int i, j;
4378	LocalHom *tmpptr;
4379	for( i=0; i<n1; i++ ) for( j=0; j<n2; j++ )
4380	{
4381	tmpptr = localhom[i][j];
4382	// fprintf( stdout, "%d-%d\n", i, j );
4383	do
4384	{
4385	fprintf( stdout, "%d-%d, reg1=%d-%d, reg2=%d-%d, imp=%f, opt=%f, wimp=%f\n", i, j, tmpptr->start1, tmpptr->end1, tmpptr->start2, tmpptr->end2, tmpptr->importance, tmpptr->opt, tmpptr->wimportance );
4386	}
4387	while( (tmpptr=tmpptr->next) );
4388	}
4389	}
4390
4391	void FreeLocalHomTable_part( LocalHom **localhomtable, int n, int m )
4392	{
4393	int i, j;
4394	LocalHom ppp, tmpptr;
4395	for( i=0; i<n; i++ )
4396	{
4397	for( j=0; j<m; j++ )
4398	{
4399	tmpptr=localhomtable[i]+j;
4400	ppp = tmpptr->next;
4401	for( ; tmpptr; tmpptr=ppp )
4402	{
4403	#if DEBUG
4404	fprintf( stderr, "i=%d, j=%d\n", i, j );
4405	#endif
4406	ppp = tmpptr->next;
4407	if( tmpptr!=localhomtable[i]+j )
4408	{
4409	#if DEBUG
4410	fprintf( stderr, "freeing %p\n", tmpptr );
4411	#endif
4412	free( tmpptr );
4413	}
4414	}
4415	}
4416	#if DEBUG
4417	fprintf( stderr, "freeing localhomtable[%d]\n", i );
4418	#endif
4419	free( localhomtable[i] );
4420	}
4421	#if DEBUG
4422	fprintf( stderr, "freeing localhomtable\n" );
4423	#endif
4424	free( localhomtable );
4425	#if DEBUG
4426	fprintf( stderr, "freed\n" );
4427	#endif
4428	}
4429
4430	void FreeLocalHomTable_two( LocalHom **localhomtable, int n, int m )
4431	{
4432	int i, j;
4433	LocalHom ppp, tmpptr;
4434	for( i=0; i<n; i++ )
4435	{
4436	for( j=0; j<m; j++ )
4437	{
4438	tmpptr=localhomtable[i]+j;
4439	ppp = tmpptr->next;
4440	for( ; tmpptr; tmpptr=ppp )
4441	{
4442	#if DEBUG
4443	fprintf( stderr, "i=%d, j=%d\n", i, j );
4444	#endif
4445	ppp = tmpptr->next;
4446	if( tmpptr!=localhomtable[i]+j )
4447	{
4448	#if DEBUG
4449	fprintf( stderr, "freeing %p\n", tmpptr );
4450	#endif
4451	free( tmpptr );
4452	}
4453	}
4454	}
4455	#if DEBUG
4456	fprintf( stderr, "freeing localhomtable[%d]\n", i );
4457	#endif
4458	free( localhomtable[i] );
4459	}
4460
4461	for( i=n; i<n+m; i++ )
4462	{
4463	for( j=0; j<n; j++ )
4464	{
4465	tmpptr=localhomtable[i]+j;
4466	ppp = tmpptr->next;
4467	for( ; tmpptr; tmpptr=ppp )
4468	{
4469	#if DEBUG
4470	fprintf( stderr, "i=%d, j=%d\n", i, j );
4471	#endif
4472	ppp = tmpptr->next;
4473	if( tmpptr!=localhomtable[i]+j )
4474	{
4475	#if DEBUG
4476	fprintf( stderr, "freeing %p\n", tmpptr );
4477	#endif
4478	free( tmpptr );
4479	}
4480	}
4481	}
4482	#if DEBUG
4483	fprintf( stderr, "freeing localhomtable[%d]\n", i );
4484	#endif
4485	free( localhomtable[i] );
4486	}
4487	#if DEBUG
4488	fprintf( stderr, "freeing localhomtable\n" );
4489	#endif
4490	free( localhomtable );
4491	#if DEBUG
4492	fprintf( stderr, "freed\n" );
4493	#endif
4494	}
4495
4496	void FreeLocalHomTable_one( LocalHom **localhomtable, int n, int m )
4497	{
4498	int i, j;
4499	LocalHom ppp, tmpptr;
4500	for( i=0; i<n; i++ )
4501	{
4502	for( j=0; j<m; j++ )
4503	{
4504	tmpptr=localhomtable[i]+j;
4505	ppp = tmpptr->next;
4506	for( ; tmpptr; tmpptr=ppp )
4507	{
4508	#if DEBUG
4509	fprintf( stderr, "i=%d, j=%d\n", i, j );
4510	#endif
4511	ppp = tmpptr->next;
4512	if( tmpptr!=localhomtable[i]+j )
4513	{
4514	#if DEBUG
4515	fprintf( stderr, "freeing %p\n", tmpptr );
4516	#endif
4517	free( tmpptr );
4518	}
4519	}
4520	}
4521	#if DEBUG
4522	fprintf( stderr, "freeing localhomtable[%d]\n", i );
4523	#endif
4524	free( localhomtable[i] );
4525	}
4526
4527	#if DEBUG
4528	fprintf( stderr, "freeing localhomtable\n" );
4529	#endif
4530	free( localhomtable );
4531	#if DEBUG
4532	fprintf( stderr, "freed\n" );
4533	#endif
4534	}
4535
4536	void FreeLocalHomTable( LocalHom **localhomtable, int n )
4537	{
4538	int i, j;
4539	LocalHom ppp, tmpptr;
4540	for( i=0; i<n; i++ )
4541	{
4542	for( j=0; j<n; j++ )
4543	{
4544	tmpptr=localhomtable[i]+j;
4545	ppp = tmpptr->next;
4546	for( ; tmpptr; tmpptr=ppp )
4547	{
4548	#if DEBUG
4549	fprintf( stderr, "i=%d, j=%d\n", i, j );
4550	#endif
4551	ppp = tmpptr->next;
4552	if( tmpptr!=localhomtable[i]+j )
4553	{
4554	#if DEBUG
4555	fprintf( stderr, "freeing %p\n", tmpptr );
4556	#endif
4557	free( tmpptr );
4558	}
4559	}
4560	}
4561	#if DEBUG
4562	fprintf( stderr, "freeing localhomtable[%d]\n", i );
4563	#endif
4564	free( localhomtable[i] );
4565	}
4566	#if DEBUG
4567	fprintf( stderr, "freeing localhomtable\n" );
4568	#endif
4569	free( localhomtable );
4570	#if DEBUG
4571	fprintf( stderr, "freed\n" );
4572	#endif
4573	}
4574
4575	char progName( char str )
4576	{
4577	char *value;
4578	if( ( value = strrchr( str, '/' ) ) != NULL )
4579	return( value+1 );
4580	else
4581	return( str );
4582	}
4583
4584	static void tabtospace( char *str )
4585	{
4586	char *p;
4587	// fprintf( stderr, "before = %s\n", str );
4588	while( NULL != ( p = strchr( str , '\t' ) ) )
4589	{
4590	*p = ' ';
4591	}
4592	// fprintf( stderr, "after = %s\n", str );
4593	}
4594
4595	static char extractfirstword( char str )
4596	{
4597	char *val = str;
4598
4599	tabtospace( str );
4600	while( *str )
4601	{
4602	if( val == str && *str == ' ' )
4603	{
4604	val++; str++;
4605	}
4606	else if( *str != ' ' )
4607	{
4608	str++;
4609	}
4610	else if( *str == ' ' )
4611	{
4612	*str = 0;
4613	}
4614	}
4615	return( val );
4616	}
4617
4618	void phylipout_pointer( FILE fp, int nseq, int maxlen, char seq, char name, int order, int namelen )
4619	{
4620	int pos, pos2, j;
4621	if( namelen == -1 ) namelen = 10;
4622	pos = 0;
4623
4624	fprintf( fp, " %d %d\n", nseq, maxlen );
4625
4626	while( pos < maxlen )
4627	{
4628	for( j=0; j<nseq; j++ )
4629	{
4630	if( pos == 0 )
4631	fprintf( fp, "%-.s", namelen, namelen, extractfirstword( name[order[j]]+1 ) );
4632	else
4633	fprintf( fp, "%-.s", namelen, namelen, "" );
4634
4635	pos2 = pos;
4636	while( pos2 < pos+41 && pos2 < maxlen )
4637	{
4638	fprintf( fp, " %.10s", seq[order[j]]+pos2 );
4639	pos2 += 10;
4640	}
4641	fprintf( fp, "\n" );
4642	}
4643	fprintf( fp, "\n" );
4644	pos += 50;
4645	}
4646	}
4647
4648	void clustalout_pointer( FILE fp, int nseq, int maxlen, char seq, char name, char mark, char comment, int order, int namelen )
4649	{
4650	int pos, j;
4651	if( namelen == -1 ) namelen = 15;
4652	pos = 0;
4653	if( comment == NULL )
4654	fprintf( fp, "CLUSTAL format alignment by MAFFT (v%s)\n\n", VERSION );
4655	else
4656	fprintf( fp, "CLUSTAL format alignment by MAFFT %s (v%s)\n\n", comment, VERSION );
4657
4658	while( pos < maxlen )
4659	{
4660	fprintf( fp, "\n" );
4661	for( j=0; j<nseq; j++ )
4662	{
4663	fprintf( fp, "%-.s ", namelen, namelen, extractfirstword( name[order[j]]+1 ) );
4664	fprintf( fp, "%.60s\n", seq[order[j]]+pos ); // Ä¹€µ€¬°ã€Š€È€À€á¡£
4665	}
4666	if( mark )
4667	{
4668	fprintf( fp, "%-.s ", namelen, namelen, "" );
4669	fprintf( fp, "%.60s\n", mark + pos ); // Ä¹€µ€¬°ã€Š€È€À€á¡£
4670	}
4671	pos += 60;
4672	}
4673	}
4674
4675
4676	void writeData_reorder_pointer( FILE fp, int locnjob, char name, int nlen, char *aseq, int order )
4677	{
4678	int i, j, k;
4679	int nalen;
4680
4681	for( i=0; i<locnjob; i++ )
4682	{
4683	k = order[i];
4684	#if DEBUG
4685	fprintf( stderr, "i = %d in writeData\n", i );
4686	#endif
4687	nalen = strlen( aseq[k] );
4688	fprintf( fp, ">%s\n", name[k]+1 );
4689	for( j=0; j<nalen; j=j+C )
4690	{
4691	#if 0
4692	strncpy( b, aseq[k]+j, C ); b[C] = 0;
4693	fprintf( fp, "%s\n",b );
4694	#else
4695	fprintf( fp, "%.*s\n", C, aseq[k]+j );
4696	#endif
4697	}
4698	}
4699	}
4700	void writeData_reorder( FILE fp, int locnjob, char name[][B], int nlen[], char aseq, int order )
4701	{
4702	int i, j, k;
4703	int nalen;
4704
4705	for( i=0; i<locnjob; i++ )
4706	{
4707	k = order[i];
4708	#if DEBUG
4709	fprintf( stderr, "i = %d in writeData\n", i );
4710	#endif
4711	nalen = strlen( aseq[k] );
4712	fprintf( fp, ">%s\n", name[k]+1 );
4713	for( j=0; j<nalen; j=j+C )
4714	{
4715	#if 0
4716	strncpy( b, aseq[k]+j, C ); b[C] = 0;
4717	fprintf( fp, "%s\n",b );
4718	#else
4719	fprintf( fp, "%.*s\n", C, aseq[k]+j );
4720	#endif
4721	}
4722	}
4723	}
4724	static void showaamtxexample()
4725	{
4726	fprintf( stderr, "Format error in aa matrix\n" );
4727	fprintf( stderr, "# Example:\n" );
4728	fprintf( stderr, "# comment\n" );
4729	fprintf( stderr, " A R N D C Q E G H I L K M F P S T W Y V\n" );
4730	fprintf( stderr, "A 4 -1 -2 -2 0 -1 -1 0 -2 -1 -1 -1 -1 -2 -1 1 0 -3 -2 0\n" );
4731	fprintf( stderr, "R -1 5 0 -2 -3 1 0 -2 0 -3 -2 2 -1 -3 -2 -1 -1 -3 -2 -3\n" );
4732	fprintf( stderr, "...\n" );
4733	fprintf( stderr, "V 0 -3 -3 -3 -1 -2 -2 -3 -3 3 1 -2 1 -1 -2 -2 0 -3 -1 4\n" );
4734	fprintf( stderr, "frequency 0.07 0.05 0.04 0.05 0.02 .. \n" );
4735	fprintf( stderr, "# Example end\n" );
4736	fprintf( stderr, "Only the lower half is loaded\n" );
4737	fprintf( stderr, "The last line (frequency) is optional.\n" );
4738	exit( 1 );
4739	}
4740
4741	double *loadaamtx( void )
4742	{
4743	int i, j, k, ii, jj;
4744	double *val;
4745	double **raw;
4746	int *map;
4747	char *aaorder = "ARNDCQEGHILKMFPSTWYV";
4748	char *inorder;
4749	char *line;
4750	char *ptr1;
4751	char *ptr2;
4752	char *mtxfname = "_aamtx";
4753	FILE *mf;
4754
4755	raw = AllocateDoubleMtx( 21, 20 );
4756	val = AllocateDoubleVec( 420 );
4757	map = AllocateIntVec( 20 );
4758
4759	if( dorp != 'p' )
4760	{
4761	fprintf( stderr, "User-defined matrix is not supported for DNA\n" );
4762	exit( 1 );
4763	}
4764
4765	mf = fopen( mtxfname, "r" );
4766	if( mf == NULL )
4767	{
4768	fprintf( stderr, "Cannot open the _aamtx file\n" );
4769	exit( 1 );
4770	}
4771
4772	inorder = calloc( 1000, sizeof( char ) );
4773	line = calloc( 1000, sizeof( char ) );
4774
4775
4776	while( !feof( mf ) )
4777	{
4778	fgets( inorder, 999, mf );
4779	if( inorder[0] != '#' ) break;
4780	}
4781	ptr1 = ptr2 = inorder;
4782	while( *ptr2 )
4783	{
4784	if( isalpha( *ptr2 ) )
4785	{
4786	ptr1 = toupper( ptr2 );
4787	ptr1++;
4788	}
4789	ptr2++;
4790	}
4791	inorder[20] = 0;
4792
4793	for( i=0; i<20; i++ )
4794	{
4795	ptr2 = strchr( inorder, aaorder[i] );
4796	if( ptr2 == NULL )
4797	{
4798	fprintf( stderr, "%c: not found in the first 20 letters.\n", aaorder[i] );
4799	showaamtxexample();
4800	}
4801	else
4802	{
4803	map[i] = ptr2 - inorder;
4804	}
4805	}
4806
4807	i = 0;
4808	while( !feof( mf ) )
4809	{
4810	fgets( line, 999, mf );
4811	// fprintf( stderr, "line = %s\n", line );
4812	if( line[0] == '#' ) continue;
4813	ptr1 = line;
4814	// fprintf( stderr, "line = %s\n", line );
4815	for( j=0; j<=i; j++ )
4816	{
4817	while( !isdigit( ptr1 ) && ptr1 != '-' && *ptr1 != '.' )
4818	ptr1++;
4819
4820	raw[i][j] = atof( ptr1 );
4821	// fprintf( stderr, "raw[][]=%f, %c-%c %d-%d\n", raw[i][j], inorder[i], inorder[j], i, j );
4822	ptr1 = strchr( ptr1, ' ' );
4823	if( ptr1 == NULL && j<i) showaamtxexample();
4824	}
4825	i++;
4826	if( i > 19 ) break;
4827	}
4828
4829	for( i=0; i<20; i++ ) raw[20][i] = -1.0;
4830	while( !feof( mf ) )
4831	{
4832	fgets( line, 999, mf );
4833	if( line[0] == 'f' )
4834	{
4835	// fprintf( stderr, "line = %s\n", line );
4836	ptr1 = line;
4837	for( j=0; j<20; j++ )
4838	{
4839	while( !isdigit( ptr1 ) && ptr1 != '-' && *ptr1 != '.' )
4840	ptr1++;
4841
4842	raw[20][j] = atof( ptr1 );
4843	// fprintf( stderr, "raw[20][]=%f, %c %d\n", raw[20][j], inorder[i], j );
4844	ptr1 = strchr( ptr1, ' ' );
4845	if( ptr1 == NULL && j<19) showaamtxexample();
4846	}
4847	break;
4848	}
4849	}
4850
4851	k = 0;
4852	for( i=0; i<20; i++ )
4853	{
4854	for( j=0; j<=i; j++ )
4855	{
4856	if( i != j )
4857	{
4858	ii = MAX( map[i], map[j] );
4859	jj = MIN( map[i], map[j] );
4860	}
4861	else ii = jj = map[i];
4862	val[k++] = raw[ii][jj];
4863	// fprintf( stderr, "%c-%c, %f\n", aaorder[i], aaorder[j], val[k-1] );
4864	}
4865	}
4866	for( i=0; i<20; i++ ) val[400+i] = raw[20][map[i]];
4867
4868	fprintf( stderr, "inorder = %s\n", inorder );
4869	fclose( mf );
4870	free( inorder );
4871	free( line );
4872	FreeDoubleMtx( raw );
4873	free( map );
4874	return( val );
4875	}
4876
4877	static void tab2space( char *s ) // nen no tame
4878	{
4879	while( *s )
4880	{
4881	if( s == '\t' ) s = ' ';
4882	s++;
4883	}
4884	}
4885
4886	static int readasubalignment( char s, int t, int *preservegaps )
4887	{
4888	int v = 0;
4889	char status = 's';
4890	char *pt = s;
4891	*preservegaps = 0;
4892	tab2space( s );
4893	while( *pt )
4894	{
4895	if( *pt == ' ' )
4896	{
4897	status = 's';
4898	}
4899	else
4900	{
4901	if( status == 's' )
4902	{
4903	if( pt == '\n' \|\| pt == '#' ) break;
4904	status = 'n';
4905	t[v] = atoi( pt );
4906	if( t[v] == 0 )
4907	{
4908	fprintf( stderr, "Format error? Sequences must be specified as 1, 2, 3...\n" );
4909	exit( 1 );
4910	}
4911	if( t[v] < 0 ) *preservegaps = 1;
4912	t[v] = abs( t[v] );
4913	t[v] -= 1;
4914	v++;
4915	}
4916	}
4917	pt++;
4918	}
4919	t[v] = -1;
4920	return( v );
4921	}
4922
4923	static int countspace( char *s )
4924	{
4925	int v = 0;
4926	char status = 's';
4927	char *pt = s;
4928	tab2space( s );
4929	while( *pt )
4930	{
4931	if( *pt == ' ' )
4932	{
4933	status = 's';
4934	}
4935	else
4936	{
4937	if( status == 's' )
4938	{
4939	if( pt == '\n' \|\| pt == '#' ) break;
4940	v++;
4941	status = 'n';
4942	if( atoi( pt ) == 0 )
4943	{
4944	fprintf( stderr, "Format error? Sequences should be specified as 1, 2, 3...\n" );
4945	exit( 1 );
4946	}
4947	}
4948	}
4949	pt++;
4950	}
4951	return( v );
4952	}
4953
4954
4955	void readsubalignmentstable( int nseq, int *table, int preservegaps, int nsubpt, int maxmempt )
4956	{
4957	FILE *fp;
4958	char *line;
4959	int linelen = 1000000;
4960	int nmem;
4961	int lpos;
4962	int i, p;
4963	int *tab01;
4964
4965	line = calloc( linelen, sizeof( char ) );
4966	fp = fopen( "_subalignmentstable", "r" );
4967	if( !fp )
4968	{
4969	fprintf( stderr, "Cannot open _subalignmentstable\n" );
4970	exit( 1 );
4971	}
4972	if( table == NULL )
4973	{
4974	*nsubpt = 0;
4975	*maxmempt = 0;
4976	while( 1 )
4977	{
4978	fgets( line, linelen-1, fp );
4979	if( feof( fp ) ) break;
4980	if( line[strlen(line)-1] != '\n' )
4981	{
4982	fprintf( stderr, "too long line? \n" );
4983	exit( 1 );
4984	}
4985	if( line[0] == '#' ) continue;
4986	if( atoi( line ) == 0 ) continue;
4987	nmem = countspace( line );
4988	if( nmem > maxmempt ) maxmempt = nmem;
4989	(*nsubpt)++;
4990	}
4991	}
4992	else
4993	{
4994	tab01 = calloc( nseq, sizeof( int ) );
4995	for( i=0; i<nseq; i++ ) tab01[i] = 0;
4996	lpos = 0;
4997	while( 1 )
4998	{
4999	fgets( line, linelen-1, fp );
5000	if( feof( fp ) ) break;
5001	if( line[strlen(line)-1] != '\n' )
5002	{
5003	fprintf( stderr, "too long line? \n" );
5004	exit( 1 );
5005	}
5006	if( line[0] == '#' ) continue;
5007	if( atoi( line ) == 0 ) continue;
5008	readasubalignment( line, table[lpos], preservegaps+lpos );
5009	for( i=0; (p=table[lpos][i])!=-1; i++ )
5010	{
5011	if( tab01[p] )
5012	{
5013	fprintf( stderr, "\nSequence %d appears in different groups.\n", p+1 );
5014	fprintf( stderr, "Hierarchical grouping is not supported.\n\n" );
5015	exit( 1 );
5016	}
5017	tab01[p] = 1;
5018	if( p > nseq-1 )
5019	{
5020	fprintf( stderr, "Sequence %d does not exist in the input sequence file.\n", p+1 );
5021	exit( 1 );
5022	}
5023	}
5024	lpos++;
5025	}
5026	free( tab01 );
5027	}
5028	fclose( fp );
5029	free( line );
5030	}
5031
5032
5033	void readmccaskill( FILE fp, RNApair *pairprob, int length )
5034	{
5035	char gett[1000];
5036	int *pairnum;
5037	int i;
5038	int left, right;
5039	float prob;
5040	int c;
5041
5042	pairnum = (int *)calloc( length, sizeof( int ) );
5043	for( i=0; i<length; i++ ) pairnum[i] = 0;
5044
5045	c = getc( fp );
5046	{
5047	if( c != '>' )
5048	{
5049	fprintf( stderr, "format error in hat4 - 1\n" );
5050	exit( 1 );
5051	}
5052	}
5053	fgets( gett, 999, fp );
5054	while( 1 )
5055	{
5056	if( feof( fp ) ) break;
5057	c = getc( fp );
5058	ungetc( c, fp );
5059	if( c == '>' \|\| c == EOF )
5060	{
5061	break;
5062	}
5063	fgets( gett, 999, fp );
5064	// fprintf( stderr, "gett = %s\n", gett );
5065	sscanf( gett, "%d %d %f", &left, &right, &prob );
5066
5067	if( left >= length \|\| right >= length )
5068	{
5069	fprintf( stderr, "format error in hat4 - 2\n" );
5070	exit( 1 );
5071	}
5072
5073	if( prob < 0.01 ) continue; // 080607, mafft ni dake eikyou
5074
5075	if( left != right && prob > 0.0 )
5076	{
5077	pairprob[left] = (RNApair )realloc( pairprob[left], (pairnum[left]+2) sizeof( RNApair ) );
5078	pairprob[left][pairnum[left]].bestscore = prob;
5079	pairprob[left][pairnum[left]].bestpos = right;
5080	pairnum[left]++;
5081	pairprob[left][pairnum[left]].bestscore = -1.0;
5082	pairprob[left][pairnum[left]].bestpos = -1;
5083	// fprintf( stderr, "%d-%d, %f\n", left, right, prob );
5084
5085	pairprob[right] = (RNApair )realloc( pairprob[right], (pairnum[right]+2) sizeof( RNApair ) );
5086	pairprob[right][pairnum[right]].bestscore = prob;
5087	pairprob[right][pairnum[right]].bestpos = left;
5088	pairnum[right]++;
5089	pairprob[right][pairnum[right]].bestscore = -1.0;
5090	pairprob[right][pairnum[right]].bestpos = -1;
5091	// fprintf( stderr, "%d-%d, %f\n", right, left, prob );
5092	}
5093	}
5094	free( pairnum );
5095	}
5096
5097	void readpairfoldalign( FILE fp, char s1, char s2, char aln1, char aln2, int q1, int q2, int of1, int *of2, int sumlen )
5098	{
5099	char gett[1000];
5100	int *maptoseq1;
5101	int *maptoseq2;
5102	char dumc;
5103	int dumi;
5104	char sinseq[100], sinaln[100];
5105	int posinseq, posinaln;
5106	int alnlen;
5107	int i;
5108	int pos1, pos2;
5109	char pa1, pa2;
5110	char qstr[1000];
5111
5112	*of1 = -1;
5113	*of2 = -1;
5114
5115	maptoseq1 = AllocateIntVec( sumlen+1 );
5116	maptoseq2 = AllocateIntVec( sumlen+1 );
5117
5118	posinaln = 0; // foldalign ga alingment wo kaesanaitok no tame.
5119
5120	while( !feof( fp ) )
5121	{
5122	fgets( gett, 999, fp );
5123	if( !strncmp( gett, "; ALIGNING", 10 ) ) break;
5124	}
5125	sprintf( qstr, "; ALIGNING %d against %d\n", q1+1, q2+1 );
5126	if( strcmp( gett, qstr ) )
5127	{
5128	fprintf( stderr, "Error in FOLDALIGN\n" );
5129	fprintf( stderr, "qstr = %s, but gett = %s\n", qstr, gett );
5130	exit( 1 );
5131	}
5132
5133	while( !feof( fp ) )
5134	{
5135	fgets( gett, 999, fp );
5136	if( !strncmp( gett, "; --------", 10 ) ) break;
5137	}
5138
5139
5140	while( !feof( fp ) )
5141	{
5142	fgets( gett, 999, fp );
5143	if( !strncmp( gett, "; ********", 10 ) ) break;
5144	// fprintf( stderr, "gett = %s\n", gett );
5145	sscanf( gett, "%c %c %s %s %d %d", &dumc, &dumc, sinseq, sinaln, &dumi, &dumi );
5146	posinaln = atoi( sinaln );
5147	posinseq = atoi( sinseq );
5148	// fprintf( stderr, "posinseq = %d\n", posinseq );
5149	// fprintf( stderr, "posinaln = %d\n", posinaln );
5150	maptoseq1[posinaln-1] = posinseq-1;
5151	}
5152	alnlen = posinaln;
5153
5154	while( !feof( fp ) )
5155	{
5156	fgets( gett, 999, fp );
5157	if( !strncmp( gett, "; --------", 10 ) ) break;
5158	}
5159
5160	while( !feof( fp ) )
5161	{
5162	fgets( gett, 999, fp );
5163	if( !strncmp( gett, "; ********", 10 ) ) break;
5164	// fprintf( stderr, "gett = %s\n", gett );
5165	sscanf( gett, "%c %c %s %s %d %d", &dumc, &dumc, sinseq, sinaln, &dumi, &dumi );
5166	posinaln = atof( sinaln );
5167	posinseq = atof( sinseq );
5168	// fprintf( stderr, "posinseq = %d\n", posinseq );
5169	// fprintf( stderr, "posinaln = %d\n", posinaln );
5170	maptoseq2[posinaln-1] = posinseq-1;
5171	}
5172	if( alnlen != posinaln )
5173	{
5174	fprintf( stderr, "Error in foldalign?\n" );
5175	exit( 1 );
5176	}
5177
5178	pa1 = aln1;
5179	pa2 = aln2;
5180	for( i=0; i<alnlen; i++ )
5181	{
5182	pos1 = maptoseq1[i];
5183	pos2 = maptoseq2[i];
5184
5185	if( pos1 > -1 )
5186	*pa1++ = s1[pos1];
5187	else
5188	*pa1++ = '-';
5189
5190	if( pos2 > -1 )
5191	*pa2++ = s2[pos2];
5192	else
5193	*pa2++ = '-';
5194	}
5195	*pa1 = 0;
5196	*pa2 = 0;
5197
5198	*of1 = 0;
5199	for( i=0; i<alnlen; i++ )
5200	{
5201	*of1 = maptoseq1[i];
5202	if( *of1 > -1 ) break;
5203	}
5204	*of2 = 0;
5205	for( i=0; i<alnlen; i++ )
5206	{
5207	*of2 = maptoseq2[i];
5208	if( *of2 > -1 ) break;
5209	}
5210
5211	// fprintf( stderr, "of1=%d, aln1 = :%s:\n", of1, aln1 );
5212	// fprintf( stderr, "of2=%d, aln2 = :%s:\n", of2, aln2 );
5213
5214	free( maptoseq1 );
5215	free( maptoseq2 );
5216	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats: