Go to Qiita Advent Calendar 2024 Top
LoginSignup
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 3 years have passed since last update.

@yoho

タンパク質の立体構造をSketchUpでモデル化 その4【主鎖を立体に】

Last updated at Posted at 2020-05-06

#はじめに
シリーズ第4回です。

このシリーズは、PDBデータをもとに、SketchUpでタンパク質の立体モデルを作成するシリーズです。
前回までに、ポリペプチドの主鎖と側鎖をedgeとして描画できるようになりました。

今回は、主鎖をfaceを持った立体(円柱が連なったもの)にして色を塗りたいと思います。コメント欄で教えていただいたカッコイイhashの使い方など(@scivola さん、ありがとうございます。)も取り入れていきたいと思います。

タンパク質の立体構造をSketchUpでモデル化 その1
タンパク質の立体構造をSketchUpでモデル化 その2【アミノ酸の側鎖を描画】
タンパク質の立体構造をSketchUpでモデル化 その3【アミノ酸の側鎖も含めて描画】
タンパク質の立体構造をSketchUpでモデル化 その4【主鎖を立体に】
タンパク質の立体構造をSketchUpでモデル化 その5【玉の描画】
タンパク質の立体構造をSketchUpでモデル化 その6【ついでにDNAをモデル化】
タンパク質の立体構造をSketchUpでモデル化 その7【最終回】

#方針
エッジの開始端にcircleを描いて、これに対してfollowmeツールを使います。この時、描画済みのedgeを指定します。

image.png

絵に描けばこんなです。circleを書く時は、中心座標と半径に加えて、circleの3次元上の向き(つまり3次元ベクトル)を指定します。向きは要するに、円を含む平面の垂線ベクトルを指定するとのことです。ここでは、1つ目のpointから2つ目のpointへと進むベクトルがちょうど垂線ベクトルになります。

Vector3dオブジェクトを得る

Point3dオブジェクトとPoint3dオブジェクトの間で引き算をすると、ベクターが得られます。

p1 = Geom::Point3d.new(1,2,3)
p2 = Geom::Point3d.new(2,2,5)

vector = p1-p2
p vector
# Vector3d(-1, 0, -2)

動作確認

円を描いて、faceにし、さらにフォローミーツールを適用します。

#モデルを取得
model = Sketchup.active_model
entities = model.active_entities

#エッジを追加
p1 = Geom::Point3d.new(1,2,3)
p2 = Geom::Point3d.new(2,2,5)
p3 = Geom::Point3d.new(5,5,5)
edges = entities.add_edges([p1,p2,p3]) #followme用に取っておく

#垂線ベクトル
vector = p1-p2

#circleを追加する。半径は適当に、0.5とする。
r = 0.5
circleEdges = entities.add_circle(p1, vector , r)

#円をfaceに変換。
face = entities.add_face(circleEdges)

#faceにフォローミーツールを適用
face.followme(edges)

動作確認です。
image.png

ちゃんと描けました。なお、垂線ベクトルをp2-p1にすると、Faceの裏地が表になった立体になってしまいます。

##色を塗る
これが難航しました。faceが1つでもあれば、それにつながっているface(間接的につながっているものも含めて)を一括して得ることのできるメソッド(all_connected)があるのですが、followmeツールを適用したfaceは、適用時に消されてしまうらしく、followmeで作った立体に対する参照を得られないのです(followmeツールの返り値はbool値)

他に良い方法があるのかもしれませんが、entitiesから一つ一つentityを取り出して、それがFaceであるかをtypenameで調べ、そうだったら先に作った円の中心座標を含むかどうかをclassify_pointで調べ、そうであったらつながっているFaceをall_connectedで取得して、それぞれに色を塗ることにしました。色を塗るには、Faceのmaterialプロパティに値("red"など)をセットします。

for x in entities
    if x.typename == "Face"
        result = x.classify_point(p1)
        if result == Sketchup::Face::PointInside
           connected_entities = x.all_connected
           for y in connected_entities
                if y.typename == "Face"
                   y.material= "red"
                end
            end
        end
    end
end

やっと目処が立ちました。

##関数の定義
Point3dオブジェクトの配列と、円柱の半径、外側の色の3つが与えられれば、立体モデルが作れることになりますので、関数 addStrand()を定義します。

# 立体的なストランドを追加するメソッド、addStrand。
#第1引数は、Point3Dの配列。ストランドの半径と色を指定
#第2引数は、円柱の半径
#第4引数は、色

def addStrand(strand,r,color)
    model = Sketchup.active_model
    entities = model.active_entities
    
    #エッジを加える
    strandEdge = entities.add_edges(strand)

    #circleを追加する。
    vector = strand[0] - strand[1]
    circleEdges = entities.add_circle(strand[0], vector , r)

    #円をfaceに変換。
    face = entities.add_face(circleEdges)
    #faceにフォローミーツールを適用
    face.followme(strandEdge)

    #ストランドに色を付ける。entitiesの中から、フォローミーに使ったFaceの中心座標を含むFaceを探す。
    for x in entities
       if x.typename == "Face"
           result = x.classify_point(strand[0])
           if result == Sketchup::Face::PointInside
               connected_entities = x.all_connected
               for y in connected_entities
                   if y.typename == "Face"
                       y.material = color
                   end
               end
           end
       end
    end
end

#コード全体
前回までのコードを多少変更して、以下のようにしました。色については、あらかじめいつくかの色の名前を配列にして保持しておき、アクセス用のindexを1ずつ増やしながらアクセスすることにしました。循環するように、indexを色の数で割り算したあまりをindexとして使います。

#モデルを取得
model = Sketchup.active_model
entities = model.active_entities

#ファイルを読み込む
lines = IO.readlines("/Users/yoho/Downloads/mmdb_5EPW.pdb")

#関数
# 立体的なストランドを追加するメソッド。ストランドの半径と色を指定
def addStrand(strand,r,color)
    model = Sketchup.active_model
    entities = model.active_entities
    
    #エッジを加える
    strandEdge = entities.add_edges(strand)
    
    #circleを追加する。
    vector = strand[0] - strand[1]
    circleEdges = entities.add_circle(strand[0], vector , r)
    
    #円をfaceに変換。
    face = entities.add_face(circleEdges)
    
    #faceにフォローミーツールを適用
    face.followme(strandEdge)
 
    #ストランドに色を付ける。entitiesの中から、フォローミーに使ったFaceの中心座標を含むFaceを探す。
    for x in entities
       if x.typename == "Face"
           result = x.classify_point(strand[0])
           if result == Sketchup::Face::PointInside
               connected_entities = x.all_connected
               for y in connected_entities
                   if y.typename == "Face"
                       y.material = color
                   end
               end
           end
       end
    end
end

#描画用情報
aminoacidsAndBonds = {
    "ALA" => [["N","CA","C"],["C","O"],["CA","CB"]],#A
    "CYS" => [["N","CA","C"],["C","O"],["CA","CB","SG"]],#C
    "ASP" => [["N","CA","C"],["C","O"],["CA","CB","CG"],["CG","OD1"],["CG","OD2"]],#D
    "GLU" => [["N","CA","C"],["C","O"],["CA","CB","CG","CD"],["CD","OE1"],["CD","OE2"]],#E
    "PHE" => [["N","CA","C"],["C","O"],["CA","CB","CG"],["CG","CD1","CE1","CZ","CE2","CD2","CG"]],#F
    "GLY" => [["N","CA","C"],["C","O"] ],#G
    "HIS" => [["N","CA","C"],["C","O"],["CA","CB","CG","ND1","CE1","NE2","CD2","CG"]],#H
    "ILE" => [["N","CA","C"],["C","O"],["CA","CB"],["CB","CG1","CD1"],["CB","CG2"]],#I
    "LYS" => [["N","CA","C"],["C","O"],["CA","CB","CG","CD","CE","NZ"]],#K
    "LEU" => [["N","CA","C"],["C","O"],["CA","CB","CG"],["CG","CD1"],["CG","CD2"]],#L
    "MET" => [["N","CA","C"],["C","O"],["CA","CB","CG","SD","CE"]],#M
    "ASN" => [["N","CA","C"],["C","O"],["CA","CB","CG"],["CG","OD1"],["CG","ND2"]],#N
    "PRO" => [["N","CA","C"],["C","O"],["CA","CB","CG","CD","N"]],#P
    "GLN" => [["N","CA","C"],["C","O"],["CA","CB","CG","CD"],["CD","OE1"],["CD","NE2"]],#Q
    "ARG" => [["N","CA","C"],["C","O"],["CA","CB","CG","CD","NE","CZ"],["CZ","NH1"],["CZ","NH2"]],#R
    "SER" => [["N","CA","C"],["C","O"],["CA","CB","OG"]],#S
    "THR" => [["N","CA","C"],["C","O"],["CA","CB"],["CB","OG1"],["CB","CG2"]],#T
    "VAL" => [["N","CA","C"],["C","O"],["CA","CB"],["CB","CG1"],["CB","CG2"]],#V
    "TRP" => [["N","CA","C"],["C","O"],["CA","CB","CG"],["CG","CD2","CE2","CZ2","CH2","CZ3","CE3","CD2"],["CG","CD1","NE1","CE2"]],#W
    "TYR" => [["N","CA","C"],["C","O"],["CA","CB","CG"],["CG","CD1","CE1","CZ","CE2","CD2","CG"],["CZ","OH"]]#Y
}
#使用する色のリスト
colors = ["SpringGreen","SlateBlue","MediumVioletRed","MediumBlue","Lime","LightSeaGreen","Indigo",
"Green","Gold","DodgerBlue","DarkOrchid","Crimson","Coral","RoyalBlue"]


#アミノ酸の主鎖データをキープ
strands = []
strandIdentity = ""
for line in lines
    dataInLine = line.split
    if dataInLine[0] == "ATOM"
        if strandIdentity != dataInLine[4]
            points = []
            strands.push(points)
            strandIdentity = dataInLine[4]
        end
        #アミノ酸バックボーンを集める。
        if ["C" ,"CA", "N"].include?(dataInLine[2])
            points.push(Geom::Point3d.new(dataInLine[6].to_f, dataInLine[7].to_f, dataInLine[8].to_f))
        end
    end
end

#主鎖を追加
colorIndex = 0
for strand in strands
    if strand.size >=2
        addStrand(strand,0.1,colors[colorIndex%14])#立体的にしないなら entities.add_edges(strand)
        colorIndex += 1
    end
end

#モノマー(アミノ酸)ごとにデータをキープ
# hashを使用。アミノ酸_strand_アミノ酸番号をkey、配列をvalueにする。配列は行全体を要素とする。
monomersAndLines = Hash.new{ |hash, key| hash[key] = [] }# 存在しないkeyを指定した時に、空の配列が返るカッコイイhash
for line in lines
    dataInLine = line.split
    if dataInLine[0] == "ATOM"
        key = dataInLine[3..5].join("_")
        monomersAndLines[key].push(line)
    end
end
    
# hashのkeyごとに、1つのアミノ酸に関する情報が含まれている。
# アミノ酸の側鎖を描画
monomersAndLines.each{|key,value|
    #最初の行を分割したアレイのindex 3からアミノ酸名を取得
    aa = value[0].split[3]
    #それぞれの行にあるATOMの名前と、座標をHashにしまう。
    atomsAndCoordinates = Hash.new
    for x in value
        data = x.split
        atomsAndCoordinates[data[2]] = Geom::Point3d.new(data[6].to_f, data[7].to_f, data[8].to_f)
    end
    #アミノ酸ごとに、描画すべき結合を、edgeとして追加する。
    #アミノ酸の原子情報が欠けている場合があるので、bondsから2つずつ取り出して、対応する座標データが両方得られる場合に限って、edgeを追加する。
    list = aminoacidsAndBonds[aa]
    slicedlist = list.slice(1,list.size-1)# 先頭の["N","CA","C"]は描画しない
    for bonds in slicedlist
        for index in 0..bonds.size - 2
            atomA = atomsAndCoordinates[bonds[index]]
            atomB = atomsAndCoordinates[bonds[index+1]]
            if  atomA != nil and atomB != nil then
                entities.add_edges([atomA,atomB])
            else
                puts "結合情報がありません" + key +" "+ aa + " "+ bonds[index] + "-" + bonds[index + 1]
            end
        end
    end
}

#描いた絵
以下のようになりました。当然、回したりできます。
image.png

#終わりに
あとは球を追加できるようにすること、ですね。
DNAについても立体構造データがPDBにあるのでそちらもついでに描画する予定です。

0
0
1

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?