Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

io.c

Visit:

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

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

Context Navigation

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

Download in other formats: