我是一名生物专业的学生,我正在尝试学习 perl、python 和 C,并在我的工作中使用这些脚本。所以,我有一个文件如下:
>sequence1
ATCGATCGATCG
>sequence2
AAAATTTT
>sequence3
CCCCGGGG
输出应如下所示,即每个序列的名称和每行中的字符数,并在文件末尾打印序列总数。
sequence1 12
sequence2 8
sequence3 8
Total number of sequences = 3
我可以让 perl 和 python 脚本工作,这是 python 脚本的例子:
#!/usr/bin/python
import sys
my_file = open(sys.argv[1]) #open the file
my_output = open(sys.argv[2], "w") #open output file
total_sequence_counts = 0
for line in my_file:
if line.startswith(">"):
sequence_name = line.rstrip('\n').replace(">","")
total_sequence_counts += 1
continue
dna_length = len(line.rstrip('\n'))
my_output.write(sequence_name + " " + str(dna_length) + '\n')
my_output.write("Total number of sequences = " + str(total_sequence_counts) + '\n')
现在,我想用 C 编写相同的脚本,这是我目前所取得的成就:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
int main(int argc, char *argv[])
{
input = FILE *fopen(const char *filename, "r");
output = FILE *fopen(const char *filename, "w");
double total_sequence_counts = 0;
char sequence_name[];
char line [4095]; // set a temporary line length
char buffer = (char *) malloc (sizeof(line) +1); // allocate some memory
while (fgets(line, sizeof(line), filename) != NULL) { // read until new line character is not found in line
buffer = realloc(*buffer, strlen(line) + strlen(buffer) + 1); // realloc buffer to adjust buffer size
if (buffer == NULL) { // print error message if memory allocation fails
printf("\n Memory error");
return 0;
}
if (line[0] == ">") {
sequence_name = strcpy(sequence_name, &line[1]);
total_sequence_counts += 1
}
else {
double length = strlen(line);
fprintf(output, "%s \t %ld", sequence_name, length);
}
fprintf(output, "%s \t %ld", "Total number of sequences = ", total_sequence_counts);
}
int fclose(FILE *input); // when you are done working with a file, you should close it using this function.
return 0;
int fclose(FILE *output);
return 0;
}
但是这段代码,当然是错误百出,我的问题是,尽管学习了很多,我仍然不能正确理解和使用内存分配和指针,所以我知道我在那部分特别有错误。如果您能对我的代码发表评论并了解它如何变成实际可用的脚本,那就太好了。顺便说一下,在我的实际数据中,每行的长度没有定义,所以我需要为此目的使用 malloc 和 realloc。
最佳答案
对于像这样的简单程序,您一次只查看短行,您不必担心动态内存分配。使用合理大小的本地缓冲区可能就足够了。
另一件事是,C 并不是特别适合快速而肮脏的字符串处理。例如,标准库中没有 strstrip
函数。您通常最终会自己实现此类行为。
一个示例实现如下所示:
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <ctype.h>
#define MAXLEN 80 /* Maximum line length, including null terminator */
int main(int argc, char *argv[])
{
FILE *in;
FILE *out;
char line[MAXLEN]; /* Current line buffer */
char ref[MAXLEN] = ""; /* Sequence reference buffer */
int nseq = 0; /* Sequence counter */
if (argc != 3) {
fprintf(stderr, "Usage: %s infile outfile\n", argv[0]);
exit(1);
}
in = fopen(argv[1], "r");
if (in == NULL) {
fprintf(stderr, "Couldn't open %s.\n", argv[1]);
exit(1);
}
out = fopen(argv[2], "w");
if (in == NULL) {
fprintf(stderr, "Couldn't open %s for writing.\n", argv[2]);
exit(1);
}
while (fgets(line, sizeof(line), in)) {
int len = strlen(line);
/* Strip whitespace from end */
while (len > 0 && isspace(line[len - 1])) len--;
line[len] = '\0';
if (line[0] == '>') {
/* First char is '>': copy from second char in line */
strcpy(ref, line + 1);
} else {
/* Other lines are sequences */
fprintf(out, "%s: %d\n", ref, len);
nseq++;
}
}
fprintf(out, "Total number of sequences. %d\n", nseq);
fclose(in);
fclose(out);
return 0;
}
很多代码都是关于强制参数和打开和关闭文件的。 (如果您将 stdin
和 stdout
与文件重定向一起使用,您可以减少很多代码。)
核心是大while
循环。注意事项:
fgets
在出错或到达文件末尾时返回NULL
。- 第一行确定行的长度,然后从末尾删除空白。
- 减少长度是不够的,最后剥离的字符串必须以空字符终止
'\0'
- 当您检查行中的第一个字符时,您应该检查一个字符,而不是一个字符串。在 C 中,单引号和双引号不能互换。
">"
是两个字符的字符串文字,'>'
和终止符'\0'
。 - 在处理字符串中的字符等可数实体时,请使用整数类型,而不是 float 。 (我在这里使用了 (signed)
int
,但是因为一行中不能有负数的字符,所以使用 unsigned 类型可能更好。) - 符号
line + 1
等同于&line[1]
。 - 我展示的代码没有检查每个序列是否总是有一个引用。我将把它作为练习留给读者。
对于初学者来说,需要跟踪的内容可能很多。对于像您这样的小型文本处理任务,Python 和 Perl 绝对更适合。
编辑:上面的解决方案不适用于长序列;它仅限于 MAXLEN
个字符。但是,如果您只需要长度而不需要序列的内容,则不需要动态分配。
这是一个更新版本,它不读取行,而是读取字符。在 '>'
上下文中,它存储了引用。否则它只是保持计数:
#include <stdlib.h>
#include <stdio.h>
#include <ctype.h> /* for isspace() */
#define MAXLEN 80 /* Maximum line length, including null terminator */
int main(int argc, char *argv[])
{
FILE *in;
FILE *out;
int nseq = 0; /* Sequence counter */
char ref[MAXLEN]; /* Reference name */
in = fopen(argv[1], "r");
out = fopen(argv[2], "w");
/* Snip: Argument and file checking as above */
while (1) {
int c = getc(in);
if (c == EOF) break;
if (c == '>') {
int n = 0;
c = fgetc(in);
while (c != EOF && c != '\n') {
if (n < sizeof(ref) - 1) ref[n++] = c;
c = fgetc(in);
}
ref[n] = '\0';
} else {
int len = 0;
int n = 0;
while (c != EOF && c != '\n') {
n++;
if (!isspace(c)) len = n;
c = fgetc(in);
}
fprintf(out, "%s: %d\n", ref, len);
nseq++;
}
}
fprintf(out, "Total number of sequences. %d\n", nseq);
fclose(in);
fclose(out);
return 0;
}
注意事项:
fgetc
从文件中读取单个字节并在文件结束时返回该字节或EOF
。在此实现中,这是唯一使用的阅读功能。- 存储引用字符串也是通过
fgetc
实现的。您也可以在跳过初始尖括号后使用fgets
。 - 计数只是读取字节而不存储它们。
n
是总计数,len
是到最后一个非空格的计数。 (您的行可能只包含没有任何尾随空格的 ACGT,因此您可以跳过空格测试并使用n
而不是len
。)
关于c - 用C获取文件中每一行的长度并写入输出文件,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/26885140/