2011-11-28 54 views
2

我有兩組PDB文件(這是不能修改的標準格式)。第一組是這樣的:如果在同一行中匹配字符串,則會對文本塊進行基於位置的替換

ATOM  18 C33 Q58 d 91  -25.677 3.886 -30.044 1.00 0.00   C 
ATOM  19 C34 Q58 d 91  -24.704 4.881 -29.447 1.00 0.00   C 
ATOM  20 C35 Q58 d 91  -23.382 4.873 -30.182 1.00 0.00   C 
ATOM  21 C8 Q58 d 91  -20.295 11.484 -33.616 1.00 0.00   C 
ATOM  22 C7 Q58 d 91  -19.198 12.305 -33.381 1.00 0.00   C 
ATOM  23 C3 Q58 d 91  -18.213 12.498 -34.383 1.00 0.00   C 

,第二個雲:

HETATM 2686 C7 589 A 1  -19.344 12.177 -33.319 1.00 25.88   C 
HETATM 2687 C8 589 A 1  -20.388 11.319 -33.511 1.00 26.31   C 
HETATM 2688 C9 589 A 1  -20.364 10.691 -34.747 1.00 26.14   C 
HETATM 2689 C10 589 A 1  -19.402 10.845 -35.729 1.00 26.34   C 
HETATM 2690 N11 589 A 1  -21.334 11.123 -32.604 1.00 26.22   N 
HETATM 2691 C12 589 A 1  -21.713 9.967 -32.081 1.00 25.65   C 

每一列是由空格可變數目,使得它的內容佔用一個特定的位置範圍分離。

第7-9列表示笛卡爾空間中的x,y,z座標。我想用所有第3列(原子類型)匹配的文件1中的座標替換文件2的座標。

例如,在該示例中,輸出文件2將是:

HETATM 2686 C7 589 A 1  -19.198 12.305 -33.381 1.00 25.88   C 
HETATM 2687 C8 589 A 1  -20.295 11.484 -33.616 1.00 26.31   C 
HETATM 2688 C9 589 A 1  -20.364 10.691 -34.747 1.00 26.14   C 
HETATM 2689 C10 589 A 1  -19.402 10.845 -35.729 1.00 26.34   C 
HETATM 2690 N11 589 A 1  -21.334 11.123 -32.604 1.00 26.22   N 
HETATM 2691 C12 589 A 1  -21.713 9.967 -32.081 1.00 25.65   C 

請注意如何座標已經改變爲前兩行(原子C7和C8)。

我嘗試過awk,但它似乎太分隔符相關,這在這個例子中不好。列3(原子類型)始終位於位置14-16,而3個座標列跨越32至54.

注意:在某些情況下,某些列可能會合並。例如,在這個例子中列5,6和合並(這也可能發生與1列2):

HETATM 2804 PG ANP A1001  23.808 17.953 28.350 1.00 52.23   P 

我的解決方案這一步(慢,但工程):

while read line ; do 
atom=$(echo "$line" | cut -c13-16) 
coord=$(grep -i "$atom" ${ligand}_${chain}_dock.tmp | cut -c32-54) 
echo "$line" | sed -r "s/^(.{31})(.{23})/\1${coord}/" >> ${ligand}_${chain}_dock.pdb 
done < ${ligand}_${chain}_ref.pdb 
+0

當字段3不匹配時,你想要做什麼? – Sorpigal

+0

什麼都沒有。只要保持原來的座標。但是,應該始終有一場比賽。請記住,第1列可以合併到第2列,第4列可以合併到第5列的某些文件中。 – mirix

+0

每個字段的確切寬度將會有所幫助。 – Sorpigal

回答

3

我大概選擇一個愚蠢的方式來解決它:玩printf語句。然而它適用於你的例子。

命令:

awk -F' *' 'NR==FNR{a[$3]=$7;b[$3]=$8;c[$3]=$9;next;}\ 
{if($3 in a)printf "%s %s %-3s %s %s %3s %11s %7s %7s %5s %s %11s\n",\ 
     $1,$2,$3,$4,$5,$6,a[$3],b[$3],c[$3],$10,$11,$12; else print $0}' file1 file2 

測試你的榜樣:

kent$ awk -F' *' 'NR==FNR{a[$3]=$7;b[$3]=$8;c[$3]=$9;next;} 
{if($3 in a)printf "%s %s %-3s %s %s %3s %11s %7s %7s %5s %s %11s\n", 
     $1,$2,$3,$4,$5,$6,a[$3],b[$3],c[$3],$10,$11,$12; else print $0}' file1 file2 
HETATM 2686 C7 589 A 1  -19.198 12.305 -33.381 1.00 25.88   C 
HETATM 2687 C8 589 A 1  -20.295 11.484 -33.616 1.00 26.31   C 
HETATM 2688 C9 589 A 1  -20.364 10.691 -34.747 1.00 26.14   C 
HETATM 2689 C10 589 A 1  -19.402 10.845 -35.729 1.00 26.34   C 
HETATM 2690 N11 589 A 1  -21.334 11.123 -32.604 1.00 26.22   N 
HETATM 2691 C12 589 A 1  -21.713 9.967 -32.081 1.00 25.65   C 
+0

+1。您不需要行尾反斜槓,也不需要顯式地「打印$ 0」 - 只需「打印」即可。默認的字段分隔符是「空格」,所以你不需要使用'-F'。 –

+0

thx,我也不需要-F'*'。 :) – Kent

+0

這個解決方案有一個問題,這就是爲什麼我通常遇到awk和PDB文件的麻煩。有時列被合併。當第2列值大於99999或第4列值大於999或9999時,通常是這種情況,例如:HETATM 2804 PG ANP A1001 23.808 17.953 28.350 1.00 52.23 P – mirix

0

這工作了。我們將OFS變量設置爲「\ t」,而不是打印每一列。

更新:在OFS變量的選項卡旁邊添加了幾個空格。這給你的輸出之間有足夠的間距。

[jaypal:~/Temp] awk -v OFS="\t " 'NR==FNR{a[$3]=$7;b[$3]=$8;c[$3]=$9;next} ($3 in a) {$7=a[$3];$8=b[$3];$9=c[$3];print $0;next} {$1=$1}1' file1 file2 
HETATM  2686  C7  589  A  1  -19.198  12.305  -33.381  1.00  25.88  C 
HETATM  2687  C8  589  A  1  -20.295  11.484  -33.616  1.00  26.31  C 
HETATM  2688  C9  589  A  1  -20.364  10.691  -34.747  1.00  26.14  C 
HETATM  2689  C10  589  A  1  -19.402  10.845  -35.729  1.00  26.34  C 
HETATM  2690  N11  589  A  1  -21.334  11.123  -32.604  1.00  26.22  N 
HETATM  2691  C12  589  A  1  -21.713  9.967  -32.081  1.00  25.65  C 
1

我對正確的字段寬度進行了猜測,但是如果它們調整正確,這應該可以正常工作。

#!/usr/bin/env bash 

file1="$1" 
file2="$2" 

fw=(7 6 4 4 4 6 9 7 9 5 16 4) 

while IFS= read -r -a f2_line ; do 
    let pos=0 
    f2_fields=() 
    for width in "${fw[@]}" ; do 
      f2_fields=("${f2_fields[@]}" "${f2_line:${pos}:${width}}") 
      let pos+=width 
    done 

    printf '%s' "${f2_fields[@]:0:6}" 
    orig=1 
    while IFS= read -r -a f1_line ; do 
      let pos=0 
      f1_fields=() 
      for width in "${fw[@]}" ; do 
        f1_fields=("${f1_fields[@]}" "${f1_line:${pos}:${width}}") 
        let pos+=width 
      done 
      if [ "${f1_fields[2]}" = "${f2_fields[2]}" ] ; then 
        orig= 
        printf '%s' "${f1_fields[@]:6:3}" 
        break 
      fi 
    done < "$file1" 
    if [ ! -z "$orig" ] ; then 
      printf '%s' "${f2_fields[@]:6:3}" 
    fi 
    printf '%s' "${f2_fields[@]:9}" 
    printf '\n' 
done < "$file2" 

這當然不是很有效率。

編輯:哎呀,不得不s/5/6 /行14.現在工作。

相關問題