1

我有一个大约 200 000 行的 pdb 文本文件。每行看起来像这样：

COMPND
SOURCE    
HETATM    1  CT  100     1     -23.207  17.632  14.543
HETATM    2  CT   99     1     -22.069  18.353  15.280
HETATM    3  OH  101     1     -21.074  18.762  14.358
HETATM    4  F   103     1     -23.816  18.483  13.675
HETATM    5  F   103     1     -24.119  17.162  15.433
HETATM    6  F   103     1     -22.680  16.591  13.841
HETATM    7  HC  104     1     -21.623  17.681  16.014
HETATM    8  HC  104     1     -22.451  19.218  15.823
HETATM    9  HO  102     1     -21.040  18.108  13.673
HETATM   10  CT  100     2      -4.340 -29.478  45.144
HETATM   11  CT   99     2      -3.051 -29.846  44.395
HETATM   12  OH  101     2      -1.968 -29.072  44.880
HETATM   13  F   103     2      -4.217 -29.778  46.464
HETATM   14  F   103     2      -5.396 -30.156  44.621
HETATM   15  F   103     2      -4.551 -28.140  45.015
HETATM   16  HC  104     2      -3.178 -29.656  43.329
HETATM   17  HC  104     2      -2.829 -30.908  44.511
HETATM   18  HO  102     2      -2.315 -28.222  45.119
HETATM   19  CT  100     3     -49.455 -17.542 -31.718
HETATM   20  CT   99     3     -49.981 -18.984 -31.736
HETATM   21  OH  101     3     -48.905 -19.897 -31.607
HETATM   22  F   103     3     -48.867 -17.273 -30.521
HETATM   23  F   103     3     -50.474 -16.668 -31.929
HETATM   24  F   103     3     -48.527 -17.405 -32.704
...

我必须将 C1 的所有第一个 CT 和 C2 的第二个 CT 更改为相同的 F1、F2、F3 和 HC 到 H1、H2。

是否可以在一个小脚本中使用 awk 和 sed 更改它们？每个 C1-C2 和 F1、F2、F3 都是同一分子（三氟乙醇 - TFE）的一部分，但有许多 TFE 分子需要定义。

所以我希望它看起来像这样：

COMPND
SOURCE    
HETATM    1  C1  100     1     -23.207  17.632  14.543
HETATM    2  C2   99     1     -22.069  18.353  15.280
HETATM    3  OH  101     1     -21.074  18.762  14.358
HETATM    4  F1  103     1     -23.816  18.483  13.675
HETATM    5  F2  103     1     -24.119  17.162  15.433
HETATM    6  F3  103     1     -22.680  16.591  13.841
HETATM    7  H1  104     1     -21.623  17.681  16.014
HETATM    8  H2  104     1     -22.451  19.218  15.823
HETATM    9  HO  102     1     -21.040  18.108  13.673
HETATM   10  C1  100     2      -4.340 -29.478  45.144
HETATM   11  C2   99     2      -3.051 -29.846  44.395
HETATM   12  OH  101     2      -1.968 -29.072  44.880
HETATM   13  F1  103     2      -4.217 -29.778  46.464
HETATM   14  F2  103     2      -5.396 -30.156  44.621
HETATM   15  F3  103     2      -4.551 -28.140  45.015
HETATM   16  H1  104     2      -3.178 -29.656  43.329
HETATM   17  H2  104     2      -2.829 -30.908  44.511
HETATM   18  HO  102     2      -2.315 -28.222  45.119
HETATM   19  C1  100     3     -49.455 -17.542 -31.718
HETATM   20  C2   99     3     -49.981 -18.984 -31.736
HETATM   21  OH  101     3     -48.905 -19.897 -31.607
HETATM   22  F1  103     3     -48.867 -17.273 -30.521
HETATM   23  F2  103     3     -50.474 -16.668 -31.929
HETATM   24  F3  103     3     -48.527 -17.405 -32.704
...

谢谢

4

2 回答 2

1

你可以awk比 更容易地使用sed，尽管我毫不怀疑sed如果你真的想要它也可以这样做。

你需要：

  • 打印字段数为 1（或 2 — 小于 3）的行。
  • 当至少有 3 列时，跟踪第 3 列中的最后一个值。
  • 如果当前列是 CT、F 或 HC 之一：
    • 如果第 3 列中的最后一个值不同，则将输入第 3 列替换为第一个字母加 1；记录你输出的是 1。
    • 否则，递增计数并输出第一个字母加上计数器。
  • 否则输出该行不变。

被翻译成awk文件中的脚本的awk.script，可能是：

NF < 3 { print; next }
$3 != "CT" && $3 != "F" && $3 != "HC" { print; next }
{ if (old_col3 != $3) { counter = 0 }
  old_col3 = $3
  $3 = substr($3, 1, 1) ++counter
  print
}

而且，当它在您的数据文件（最初命名为data）上运行时，我得到：

$ awk -f awk.script data
COMPND
SOURCE    
HETATM 1 C1 100 1 -23.207 17.632 14.543
HETATM 2 C2 99 1 -22.069 18.353 15.280
HETATM    3  OH  101     1     -21.074  18.762  14.358
HETATM 4 F1 103 1 -23.816 18.483 13.675
HETATM 5 F2 103 1 -24.119 17.162 15.433
HETATM 6 F3 103 1 -22.680 16.591 13.841
HETATM 7 H1 104 1 -21.623 17.681 16.014
HETATM 8 H2 104 1 -22.451 19.218 15.823
HETATM    9  HO  102     1     -21.040  18.108  13.673
HETATM 10 C1 100 2 -4.340 -29.478 45.144
HETATM 11 C2 99 2 -3.051 -29.846 44.395
HETATM   12  OH  101     2      -1.968 -29.072  44.880
HETATM 13 F1 103 2 -4.217 -29.778 46.464
HETATM 14 F2 103 2 -5.396 -30.156 44.621
HETATM 15 F3 103 2 -4.551 -28.140 45.015
HETATM 16 H1 104 2 -3.178 -29.656 43.329
HETATM 17 H2 104 2 -2.829 -30.908 44.511
HETATM   18  HO  102     2      -2.315 -28.222  45.119
HETATM 19 C1 100 3 -49.455 -17.542 -31.718
HETATM 20 C2 99 3 -49.981 -18.984 -31.736
HETATM   21  OH  101     3     -48.905 -19.897 -31.607
HETATM 22 F1 103 3 -48.867 -17.273 -30.521
HETATM 23 F2 103 3 -50.474 -16.668 -31.929
HETATM 24 F3 103 3 -48.527 -17.405 -32.704
$

这不会保留修改后的行中的所有间距，而是您需要的。如果确实需要保留间距，则必须编写一条printf()语句来正确格式化字段（代替print最后一段代码中的 ：

printf("%s %4s %3s %4s %5s %11s %7s %7s\n", $1, $2, $3, $4, $5, $6, $7, $8);

这确实保留了间距，但通常会使代码的健壮性降低。它利用了比n in短的字符串%ns是右对齐的属性。这会产生：

COMPND
SOURCE    
HETATM    1  C1  100     1     -23.207  17.632  14.543
HETATM    2  C2   99     1     -22.069  18.353  15.280
HETATM    3  OH  101     1     -21.074  18.762  14.358
HETATM    4  F1  103     1     -23.816  18.483  13.675
HETATM    5  F2  103     1     -24.119  17.162  15.433
HETATM    6  F3  103     1     -22.680  16.591  13.841
HETATM    7  H1  104     1     -21.623  17.681  16.014
HETATM    8  H2  104     1     -22.451  19.218  15.823
HETATM    9  HO  102     1     -21.040  18.108  13.673
HETATM   10  C1  100     2      -4.340 -29.478  45.144
HETATM   11  C2   99     2      -3.051 -29.846  44.395
HETATM   12  OH  101     2      -1.968 -29.072  44.880
HETATM   13  F1  103     2      -4.217 -29.778  46.464
HETATM   14  F2  103     2      -5.396 -30.156  44.621
HETATM   15  F3  103     2      -4.551 -28.140  45.015
HETATM   16  H1  104     2      -3.178 -29.656  43.329
HETATM   17  H2  104     2      -2.829 -30.908  44.511
HETATM   18  HO  102     2      -2.315 -28.222  45.119
HETATM   19  C1  100     3     -49.455 -17.542 -31.718
HETATM   20  C2   99     3     -49.981 -18.984 -31.736
HETATM   21  OH  101     3     -48.905 -19.897 -31.607
HETATM   22  F1  103     3     -48.867 -17.273 -30.521
HETATM   23  F2  103     3     -50.474 -16.668 -31.929
HETATM   24  F3  103     3     -48.527 -17.405 -32.704

由于看起来当您达到 10,000 条记录时，HETATM列和后面的数字列合并为一列：

HETATM   21  OH  101     3     -48.905 -19.897 -31.607
…
HETATM 9999  HO  102  1111     -24.504 -16.257 -35.613
HETATM10000  CT  100  1112       9.045  23.978  29.038
HETATM10001  CT   99  1112      10.488  24.501  29.083
HETATM10002  OH  101  1112      11.370  23.545  28.522
HETATM10003  F   103  1112       8.650  23.804  27.749
HETATM10004  F   103  1112       8.209  24.855  29.654
HETATM10005  F   103  1112       8.996  22.779  29.679

目前尚不清楚当数字达到 100,000 及以上时会发生什么。但是，可以通过计算列数并适当地工作来处理它（在大多数情况下）。

NF < 7 { print; next }
NF == 8 && $3 != "CT" && $3 != "F" && $3 != "HC" { print; next }
NF == 7 && $2 != "CT" && $2 != "F" && $2 != "HC" { print; next }
NF == 8 {
          if (old_mark != $3) { counter = 0 }
          old_mark = $3
          $3 = substr($3, 1, 1) ++counter
          printf("%s %4s %3s %4s %5s %11s %7s %7s\n", $1, $2, $3, $4, $5, $6, $7, $8);
        }
NF == 7 {
          if (old_mark != $2) { counter = 0 }
          old_mark = $2
          $2 = substr($2, 1, 1) ++counter
          printf("%s %3s %4s %5s %11s %7s %7s\n", $1, $2, $3, $4, $5, $6, $7);
        }

请注意使用“列号中性”名称old_mark。如果第 9,999 行包含CT并且第 10,000 行也包含CT，则映射需要是连续的 (C1, C2) 等。您可以使用：

NF < 7 { print; next }
NF == 8 && $3 != "CT" && $3 != "F" && $3 != "HC" { print; next }
NF == 7 && $2 != "CT" && $2 != "F" && $2 != "HC" { print; next }
{
    colnum = NF - 5
    if (old_mark != $colnum) { counter = 0 }
    old_mark = $colnum
    $colnum = substr($colnum, 1, 1) ++counter
    if (NF == 7)
        printf("%s %3s %4s %5s %11s %7s %7s\n", $1, $2, $3, $4, $5, $6, $7);
    else
        printf("%s %4s %3s %4s %5s %11s %7s %7s\n", $1, $2, $3, $4, $5, $6, $7, $8);
}

可能有一种方法可以使用一个printf()电话，但我怀疑它是否值得付出努力。

于 2015-05-14T21:53:33.770 回答
0

while read这是使用循环解决此问题的一种方法，grep并且sed：

counter=0
while read line; do 
  # if a line has CT, CF, F in it... 
  if echo $line | grep -Pq '(CT|HC|F) '; then 
    # increment the counter and...
    counter=$((counter+1))

    # replace the 15th character with the counter!
    echo $line | sed "s/./$counter/15"
  else 

    # otherwise, reset the counter, and echo the line
    counter=0
    echo $line
  fi
done < molecule.txt

然后，您可以将其通过管道传输到另一个文件或标准输出！

于 2015-05-14T21:12:12.577 回答