Mini Map Aide (gu2ifv) - PasteCode.io

Mini Map Aide

File for https://minimapai.de

jbcooper

Public2 years ago4 Snippets

index.html

luaa year ago

index.html for minimapai.de

<!DOCTYPE html><html><head><meta name="viewport" content="width=device-width, initial-scale=1.0"><meta charset="UTF-8">
<title>mini map aide - display protein electron density maps on mobile phone browsers - PDB and CCP4 map file viewer for proteins, nucleotides and ligands</title>
<script src='https://minimapai.de/fengari_web.js' type="text/javascript" async></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/three.js/0.155.0/three.min.js"></script>
<script src="https://minimapai.de/surfacenets.js"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjs/10.1.0/math.min.js"></script>
<style>
  #topbar {background-color:dodgerblue; position:relative; line-height:30px; font-size:22px; text-align:center; color:white;}
  body {margin:auto; font-size:16px; font-family:arial; color:white; background-color:#0C0A0A;}
  p {margin:0px;}
  abbr{border-bottom:none !important; cursor:inherit !important; text-decoration:none !important;}
button:active {color:gray;}
</style></head><body>
<div id="topbar"><p>
<abbr title="Link to help file."><a href="minimaphelp.html" target="_blank"><button>About Mini Map Aide</button></a></abbr>
<abbr title="The start point in the sequence. This is an internal index for the residues and is not necessarily the same as the residue number in the PDB file."><label for="position" style="margin-left: 10px;">
Go&nbsp;to: </label><input type="number" min="-999" max="9999" step="1" value="1" id="position" name="position" size="4"></abbr>
<abbr title="The contour level in sigma or rms units. Change sign with +/- button."><label for="densign" style="margin-left: 10px;">
Contour&nbsp;level: <button id="densign" style="margin-left: 10px;">+/-</button></label><input type="number" min="-20" max="20" step="0.1" value="1" id="contour" name="contour" size="4"></abbr>
<abbr title="Select the PDB file."><span style="display:inline-block;"><label for="pdbInput" style="margin-left: 10px;">PDB: </label><input type="file" id="pdbInput" style="width:180px;"></span></abbr>
<abbr title="Select CCP4 map or MTZ file from refinement."><span style="display:inline-block;"><label for="mapInput" style="margin-left: 10px;">Map: </label><input type="file" id="mapInput" style="width:180px;" disabled></span></abbr>
<span style="display:inline-block;"><abbr title="Click to view outer part of sidechain or nucleotide base."><button id="viewsidechain" style="margin-left: 10px;"  disabled>Side chain</button></abbr>
<abbr title="Swap between animation and scrolling."><button id="stop" style="margin-left: 10px;"  disabled>Scroll</button></abbr>
<abbr title="Move to previous residue."><button id="back" style="margin-left: 10px;" disabled>Last residue</button></abbr>
<abbr title="Move to next residue."><button id="step" style="margin-left: 10px;" disabled>Next residue</button></abbr></span></div>
<div id="graphics" style="text-align: center">
<div><span id="info">Select a PDB file</span>
<span style="position: absolute; left: 0;"><abbr title="Tick this box to display symmetry mates. You need to click the Update button in order for this to take effect. The first time this option is used will be slow as the symmetry mates are calculated but it will be quicker after that."><input type="checkbox" id="symmetry" style="margin-top:20px; margin-left: 5px;" disabled>Symmetry</abbr>
<abbr title="Mutate."><select name="mutate" id="mutate" style="margin-top:20px; margin-left: 5px;" disabled><option value="">Mutate</option><option value="ALA">ALA</option><option value="ARG">ARG</option><option value="ASN">ASN</option><option value="ASP">ASP</option><option value="CYS">CYS</option><option value="GLN">GLN</option><option value="GLU">GLU</option><option value="GLY">GLY</option><option value="HIS">HIS</option><option value="ILE">ILE</option><option value="LEU">LEU</option><option value="LYS">LYS</option><option value="MET">MET</option><option value="PHE">PHE</option><option value="PRO">PRO</option><option value="SER">SER</option><option value="THR">THR</option><option value="TRP">TRP</option><option value="TYR">TYR</option><option value="VAL">VAL</option></select></abbr></span>
<span style="position: absolute; right: 0;"><abbr title="Click to move atoms in the residue at the centre of the display. The moveable atoms will then appear green and can be clicked and dragged. The selected atom is enlarged momentarily to indicate which one is moving."><button id="move" style="margin-top:20px; margin-right: 5px;" disabled>Move</button></abbr><abbr title="Tidy up the geometry of the residue after the atoms have been moved. This is a bit approximate and does not use the electron density map at all. It may be repeated if necessary."><button id="tidy" style="margin-right: 5px;" disabled>Tidy</button></abbr><abbr title="Undo the move and tidy operations on this residue."><button id="undo" style="margin-right: 5px;" disabled>Undo</button></abbr><abbr title="Save the changes in a new PDB file. It is *important* to do this before closing the browser tab if any of the changes made during the session are still wanted for future use!"><button id="save" style="margin-right: 5px;" disabled>Save</button></abbr></span></div>
<div id="myCanvas" style="padding: 0; max-width:100%; height:auto;"></div></div>
<script type="application/lua">

--  |d|e|f|a|u|l|t|s| |f|o|r| |m|i|n|i| |m|a|p| |a|i|d|e|                       v11feb25

mapcolpos="dodgerblue" mapcolneg="red" bondcol="#FFD633" altbondcol="#0FFF50" maxbondlen=2.0 maxssbondlen=3.0 boxrad=6.5 maprad=4.2

local js = require "js"
local window = js.global
local threejs = window.THREE
local mathjs = window.math
local document = window.document
canvas=document:getElementById("myCanvas")

el1=mathjs:index(0)  el2=mathjs:index(1)  el3=mathjs:index(2)

width=window.innerWidth
height=window.innerHeight
ratio=width/height
renderer=js.new(threejs.WebGLRenderer)
camera=js.new(threejs.PerspectiveCamera,20,ratio,30,50)
--camera.position:set(0,0,10)
camera.position.z=40 
--camp=camera.position print(camp.x,camp.y,camp.z)
scene=js.new(threejs.Scene)
scene:add(camera)
renderer:setSize(width,height)
canvas:appendChild(renderer.domElement)

material=js.new(threejs.LineBasicMaterial)
material.color:set(bondcol)
material.linewidth=2
altmat=js.new(threejs.LineBasicMaterial)
altmat.color:set(altbondcol)
altmat.linewidth=2
matermap=js.new(threejs.LineBasicMaterial)
matermap.linewidth=1

light=js.new(threejs.PointLight,0xFFFFFF)
light.position:set(50, 50, 50);
scene:add(light)
xworld = js.new(threejs.Vector3, 1, 0, 0 )
yworld = js.new(threejs.Vector3, 0, 1, 0 )
mouseVector = js.new(threejs.Vector2)
projector = js.new(threejs.Raycaster)

mouse_x,mouse_y,twang=0,0,0 moving=false spin=true viewSideCh=false foundSideCen=false sideAtGrk={"E","D","G","B"}  --ncyc=0 maxdiff=100  --last 2 normally commented out

function mouseDown(event)
MouseDown = true
mouseX = event.offsetX
mouseY = event.offsetY
get_moving_atom()
end

function touchStart(event)
event:preventDefault() -- needed to stop browser reading swipe events
MouseDown = true
mouseX = event.touches[0].pageX - event.touches[0].target.offsetLeft  -- custom event.offsetX,Y needed for mobile
mouseY = event.touches[0].pageY - event.touches[0].target.offsetTop
get_moving_atom()
end

function mouseMove(event)
if MouseDown==true then
deltaX = (event.offsetX - mouseX)
deltaY = (event.offsetY - mouseY)
if moving and intersect[0]~=nil then
	mouse_x=5*(deltaX)/width
	mouse_y=5*(-deltaY)/height
	dist=math.sqrt(mouse_x^2+mouse_y^2)
	if dist>0.1 then update_orig_coords() mouseX=event.offsetX mouseY=event.offsetY end
else
group:rotateOnWorldAxis(xworld,deltaY/20000)  --Swapped deliberately. Rotating the group on world axes gives more predictable mouse/touch behaviour.
group:rotateOnWorldAxis(yworld,deltaX/20000)
end
renderer:render(scene, camera)
end
end

function touchMove(event)
if MouseDown==true then
event:preventDefault()  -- needed to stop browser reading swipe events
deltaX = (event.touches[0].pageX - event.touches[0].target.offsetLeft - mouseX)/200  -- custom event.offsetX,Y needed for mobile
deltaY = (event.touches[0].pageY - event.touches[0].target.offsetTop - mouseY)/200
if moving and intersect[0]~=nil then
	mouse_x=500*(deltaX)/width
	mouse_y=500*(-deltaY)/height
	dist=math.sqrt(mouse_x^2+mouse_y^2)
	if dist>0.1 then update_orig_coords() mouseX=event.touches[0].pageX-event.touches[0].target.offsetLeft mouseY=event.touches[0].pageY-event.touches[0].target.offsetTop end
else
group:rotateOnWorldAxis(xworld,deltaY/40)  -- swapped deliberately
group:rotateOnWorldAxis(yworld,deltaX/40)
end
renderer:render(scene, camera)
end
end

function mouseUp(event)
MouseDown = false
if moving and intersect[0]~=nil then update_orig_coords() end
end

function touchEnd(event) mouseUp(event) end


function get_moving_atom()
mouseVector.x = 2 * (mouseX / width) - 1
mouseVector.y = 1 - 2 * ( mouseY / height )
--print(mouseVector.x,mouseVector.y)
if moving then
projector:setFromCamera( mouseVector, camera )
intersect = projector:intersectObject( moveables, true ) --only moveables are clickable
if intersect[0]~=nil then
		--print("intersect object x:",intersect[0].object.position.x)
		movatoms=intersect[0].object
		movatom=movatoms.name --print(movatom)
		--print("movatoms:",window.console:log(movatoms))
		--movatoms.material.color:set("white")
		movatoms.material.size=20
end
end
end


function startmove()
moving=true
document:getElementById("tidy").disabled=false
document:getElementById("undo").disabled=false
document:getElementById("save").disabled=false
resisafe={} for i,j in pairs(protein[startpos]) do resisafe[i]={} for k,l in pairs(j) do resisafe[i][k]=l end end	
move()
end


function move()
moveables = js.new(threejs.Group)
vecmov=js.new(threejs.Vector3,0,0,0)
geomov=js.new(threejs.BufferGeometry)
for k,l in pairs(protein[startpos]) do
	matermov=js.new(threejs.PointsMaterial)
	matermov.size=8
	matermov.sizeAttenuation=false
	matermov.color:set("#0FFF50")
	vecmove=window:Array(vecmov:set(l.x-xcent,l.y-ycent,l.z-zcent))
	geomove=geomov:clone()
	geomove:setFromPoints(vecmove)
	meshmove=js.new(threejs.Points,geomove,matermov)
	meshmove.name=k
	moveables:add(meshmove)
end
group:add(moveables)
renderer:render(scene, camera)
geomov:dispose() geomove:dispose() matermov:dispose() --RAM housekeeping
end


function update_orig_coords()
rotamat=group.matrix  --grab the group rotation matrix
rotmat=rotamat.elements   --N.B. three.js uses the transposed rotation matrix.
invrot1=window:Array(rotmat[0],rotmat[1],rotmat[2])
invrot2=window:Array(rotmat[4],rotmat[5],rotmat[6])
invrot3=window:Array(rotmat[8],rotmat[9],rotmat[10])
invrot=window:Array(invrot1,invrot2,invrot3)
transprot=js.new(mathjs.matrix,invrot)
mouarr=window:Array(mouse_x,mouse_y,0)
munrot=mathjs:multiply(transprot,mouarr)
munrotx=(mathjs:subset(munrot,el1))
munroty=(mathjs:subset(munrot,el2))
munrotz=(mathjs:subset(munrot,el3))
--print("Unrotated shifts:",munrotx,munroty,munrotz)
protein[startpos][movatom].x=protein[startpos][movatom].x+munrotx
protein[startpos][movatom].y=protein[startpos][movatom].y+munroty
protein[startpos][movatom].z=protein[startpos][movatom].z+munrotz
group:remove(moveables)
move()
end


function save()
format_pdb_normal="%s%5i  %-4s%3s%2s%-4s   %8.3f%8.3f%8.3f%s\n"
format_pdb_metals="%s%5i %-4s %3s%2s%-4s   %8.3f%8.3f%8.3f%s\n"
pdbOutputText=crystline
nAt=1
for i,res in ipairs(protein) do
for j,atom in pairs(res) do
if #j==3 then format_pdb=format_pdb_normal else format_pdb=format_pdb_metals end  --handle left shift of metals
pdbOutputLine=string.format(format_pdb,"ATOM  ",nAt,j,atom.res,atom.chain,atom.num,atom.x,atom.y,atom.z,"  1.00 10.00")
pdbOutputText=pdbOutputText..pdbOutputLine
nAt=nAt+1
end
end
a=document:createElement("a") --create imaginary download button
content=window:Array(pdbOutputText)
contentType=js.new(window.Object)
contentType["type"]="text/plain"
savefile=js.new(window.Blob,content,contentType)
a.href=window.URL:createObjectURL(savefile)
a.download=pdbinput.files[0].name
a:click()
--print("a:",window.console:log(a))
window.URL:revokeObjectURL(a.href)
moving=false
rotamat=group.matrix  --grab the group rotation matrix
scene:remove(group)
readform()
getBoxAtoms()
getMapBox()
group:add(meshmap)
renderer:render(scene, camera)
document:getElementById("save").disabled=true
end


function undo()
protein[startpos]=resisafe
run()
end


function run() --redraw
document:getElementById("info").innerHTML="Please wait"
if group~=nil then window:setTimeout(function() run2() end, 1) end --1 msec delay to let the DOM display the message
end

function run2()
rotamat=group.matrix  --grab the group rotation matrix
scene:remove(group)
readform()
getBoxAtoms()
if minimap~=nil then getMapBox() end
document:getElementById("mutate").disabled=false
document:getElementById("move").disabled=false
document:getElementById("tidy").disabled=true
document:getElementById("undo").disabled=true	
document:getElementById("save").disabled=true	
end


function rerun()  --redraw after tidy
rotamat=group.matrix  --grab the group rotation matrix
scene:remove(group)
readform()
getBoxAtoms()
if minimap~=nil then getMapBox() end
document:getElementById("save").disabled=false
end


function step(sign)
rotamat=group.matrix  --grab the group rotation matrix
scene:remove(group)
document:getElementById("tidy").disabled=true
document:getElementById("undo").disabled=true	
document:getElementById("save").disabled=true
oldstartpos=startpos
startpos=startpos+sign*1
if startpos<1 or startpos>#protein then startpos=oldstartpos end
startRes=protein[startpos]
if startRes~=nil then _,firstAtom=next(startRes)
    firstAtomResNum=string.sub(firstAtom.num,1,-2)  firstAtomResNum=tonumber(firstAtomResNum) --remove insertion code so resnum is a number
    position.value=firstAtomResNum
    chainbox.value=firstAtom.chain
end
getBoxAtoms()
if minimap~=nil then getMapBox() end
end


function stop()  --stop/start animation
if spin==true then 
spin=false
renderer:setAnimationLoop()
document:getElementById("stop").innerHTML="Roll"
else 
spin=true
renderer:setAnimationLoop(twirl)
document:getElementById("stop").innerHTML="Scroll"	
end
end


function readfile()
scene:remove(group)
pdbfile=pdbinput.files[0]
--print(pdbfile.name)
reader=js.new(window.FileReader)
--reader.onload=function() print(reader.result) end
reader.onload=readatoms
reader:readAsText(pdbfile)	
end


function readcell()
a=crystline:sub(7,15)   a=tonumber(a)
b=crystline:sub(16,24)  b=tonumber(b)
c=crystline:sub(25,33)  c=tonumber(c)
alpha_deg=crystline:sub(34,40)  alpha_deg=tonumber(alpha_deg)
beta_deg=crystline:sub(41,47)  beta_deg=tonumber(beta_deg)
gamma_deg=crystline:sub(48,54)  gamma_deg=tonumber(gamma_deg)
if a==nil or b==nil or c==nil or alpha_deg==nil or beta_deg==nil or gamma_deg==nil then print("Incomplete unit cell parameters.") readcell_ok=false --^&
else
    sg=crystline:sub(56,65)
    print("Space group: "..sg.."Unit cell: a="..a..", b="..b..", c="..c..", alpha="..alpha_deg..", beta="..beta_deg..", gamma="..gamma_deg)
    if sg=="P 1       " then 
		symmetry_ok()
    else
        symmetry=symops[sg]
        if symmetry==nil then print("Space group not recognised.") readcell_ok=false
        else 
		sg_rotations=symmetry.rotas
		sg_translations=symmetry.trans		
		symmetry_ok() end
    end
end  --^&
end


function symmetry_ok() fractorth()  readcell_ok=true  document:getElementById("symmetry").disabled=false  end


function fractorth()
alpha=math.rad(alpha_deg) beta=math.rad(beta_deg) gamma=math.rad(gamma_deg)
cos_alpha=math.cos(alpha)  cos_beta=math.cos(beta)  cos_gamma=math.cos(gamma)
cos_alpha_star=(cos_beta*cos_gamma-cos_alpha)/(math.sin(beta)*math.sin(gamma))
sin_alpha_star=math.sqrt(1-cos_alpha_star^2)
cos_delta=-math.sin(beta)*cos_alpha_star
cos_epsilon=math.sin(beta)*sin_alpha_star
L1=window:Array(a,b*cos_gamma,c*cos_beta)
L2=window:Array(0,b*math.sin(gamma),-c*math.sin(beta)*cos_alpha_star)
L3=window:Array(0,0,c*math.sin(beta)*sin_alpha_star)
La=window:Array(L1, L2, L3)
L=js.new(mathjs.matrix,La)
frorthlog="Orthogonalisation matrix\n"
for i=0,2 do for j=0,2 do frorthlog=frorthlog.." "..string.format("%14.9f",L:subset(mathjs:index(i, j))) end frorthlog=frorthlog.."\n" end
Linv=mathjs:inv(L)  --fractionalisation matrix
frorthlog=frorthlog.."Fractionalisation matrix\n"
for i=0,2 do for j=0,2 do frorthlog=frorthlog.." "..string.format("%14.9f",Linv:subset(mathjs:index(i, j))) end frorthlog=frorthlog.."\n" end
print(frorthlog)
end


function readatoms()  --reads all atoms into table protein={res1,res2,..}
nAt,nCA,oldresnum,xtot,ytot,ztot=0,0,0,0,0,0  protein,residue,chains={},{},{}  chainold=""  --residue={["CA "]={["res"]=resnam,["num"]=resnum,["x"]=3.1,["y"]=..},["CB "]={},..}
Xmin,Ymin,Zmin,Xmax,Ymax,Zmax=math.huge,math.huge,math.huge,-math.huge,-math.huge,-math.huge
filetext=reader.result
for line in filetext:gmatch("[^\r\n]+") do
	if line:sub(1,6)=="CRYST1" then crystline=line.."\n" readcell()
	elseif line:sub(1,6)=="ATOM  " or line:sub(1,6)=="HETATM" then
		atnam=line:sub(14,16)
		altgrp=line:sub(17,17) --print(altgrp)
		attype=atnam:sub(1,1)
		metal=line:sub(13,13) if metal~=" " then  if metal=="C" then metal="X" end  atnam=metal..atnam attype=metal..attype end  --remove left shift of heavier elements and long H-atom names, replace "C" in metals with "X" to allow atom drawing
			--if attype~="H" and attype~="G" and attype~="D" and attype~="E" and attype~="Z" then  --omit H's and long H-atom names
			if attype~="H" then  --omit H atoms
			resnam=line:sub(18,20)
			chain=line:sub(21,22)
			resnum=line:sub(23,27)
			if resnum~=oldresnum or chain~=chainold then oldresnum=resnum if nAt>0 then table.insert(protein,residue) residue={} end end
			if chain~=chainold then
				uniqch=true
				for i,j in ipairs(chains) do if chain==j then uniqch=false end end
				if uniqch then table.insert(chains,chain) print("Chain:", chain) end
				chainold=chain
			end
			x_i=tonumber(line:sub(28,38))
			y_i=tonumber(line:sub(39,46))
			z_i=tonumber(line:sub(47,54))
			x_o=x_i y_o=y_i z_o=z_i
			if atnam=="CA " then xtot=xtot+x_o  ytot=ytot+y_o  ztot=ztot+z_o nCA=nCA+1 end
			if x_o>Xmax then Xmax=x_o end  if y_o>Ymax then Ymax=y_o end  if z_o>Zmax then Zmax=z_o end
			if x_o<Xmin then Xmin=x_o end  if y_o<Ymin then Ymin=y_o end  if z_o<Zmin then Zmin=z_o end
			atom={["res"]=resnam,["num"]=resnum,["chain"]=chain,["x"]=x_o,["y"]=y_o,["z"]=z_o,["nsymop"]={0},["postrot_shift"]={0,0,0}}
			if altgrp~=" " then atnam=atnam..altgrp end
			residue[atnam]=atom
			nAt=nAt+1
			if nAt==1 then firstAtomResNum=resnum:sub(1,-2)  firstAtomResNum=tonumber(firstAtomResNum) end
		end
	end
end
table.insert(protein,residue) -- this line only applies to last residue
--print("Number of atoms read: ", n," Last residue number: ",resnum,"\n")
document:getElementById("stop").disabled=false	
document:getElementById("step").disabled=false
document:getElementById("back").disabled=false
document:getElementById("viewsidechain").disabled=false	
document:getElementById("mapInput").disabled=false
--document:getElementById("info").innerHTML="Select CCP4 map"
position:insertAdjacentHTML("beforebegin","<input id=\"chain\" name=\"chain\" size=\"2\">")
chainbox=document:getElementById("chain")
chainbox:addEventListener("change", run)  --add callback to the chain box
--document:getElementById("position").value=firstAtomResNum
--for i,res in ipairs(protein) do print("Residue:",i) for j,k in pairs(res) do print("Atom:",j) for l,m in pairs(k) do print(l,m) end end end
--for i,res in ipairs(protein) do _,atm=next(res) print(atm.chain,atm.num) end
getBoxAtoms()
end


function getCentre()
print("Num residues, current index:",#protein,startpos)
if protein[startpos]==nil then getStartPos() end

_,firstAtom=next(protein[startpos])
firstAtomResNum=string.sub(firstAtom.num,1,-2)  firstAtomResNum=tonumber(firstAtomResNum) --remove insertion code so resnum is a number
position.value=firstAtomResNum
chainbox.value=firstAtom.chain

if protein[startpos]["CA "]~= nil then centrat=protein[startpos]["CA "]
--else for i,j in pairs(protein[startpos]) do centrat=protein[startpos][i] break end -- if no CA, iterator takes first atom
else centrat=firstAtom  -- if no CA, next takes first atom
end
if viewSideCh==true then
for m,n in ipairs(sideAtGrk) do
    for i,j in pairs(protein[startpos]) do
        if i:sub(2,2)==n or i:sub(2,3)=="5 " then  -- go up to 4 steps along a side chain or find atom 5 of base
           centrat=protein[startpos][i]
           foundSideCen=true
           break
        end
    end
if foundSideCen==true then break end
end
end
--print(centrat.x,centrat.y,centrat.z)
xcent,ycent,zcent=centrat.x,centrat.y,centrat.z
xmin,xmax,ymin,ymax,zmin,zmax=centrat.x-boxrad,centrat.x+boxrad,centrat.y-boxrad,centrat.y+boxrad,centrat.z-boxrad,centrat.z+boxrad
document:getElementById("info").innerHTML=centrat.res.." "..centrat.num..centrat.chain
end


function getBonds(thing)
box={}
for i,res in ipairs(thing) do  --##
residue={}
for j,atom in pairs(res) do
	--print(j,atom.x)
	if xmin<=atom.x and atom.x<=xmax and ymin<=atom.y and atom.y<=ymax and zmin<=atom.z and atom.z<=zmax then
		--print("Inside box:",atom.res,atom.num,j)
		residue[j]=atom
	end
end
if next(residue)~=nil then table.insert(box,residue) end
end

extrabox={}  -- now only for Cys SG atoms
for i,res in ipairs(box) do  --##  gets the usual bonds
bonds={} cys={} jj=0
for j,atom in pairs(res) do
	jj=jj+1	kk=0 --kk must be reset here (bug found 11aug23)
	if j=="SG " then cys["SG "]=atom table.insert(extrabox,cys) end  -- SG's added to extrabox
	for k,atom2 in pairs(res) do
                other_part=false
                kk=kk+1
		if kk>jj then
                        --print("j=",j,"k=",k,"jj,kk",jj,kk)
                        --if #k==4 and #j==4 then print("kj=",k[4],j[4]) if string.sub(k,4)~=string.sub(j,4) then other_part=true end end
                        if #k==4 and #j==4 then if string.sub(k,4)~=string.sub(j,4) then other_part=true end end
                        d=math.sqrt((atom.x-atom2.x)^2+(atom.y-atom2.y)^2+(atom.z-atom2.z)^2)
			if not other_part and d<maxbondlen then
			bond={{atom.x,atom.y,atom.z},{atom2.x,atom2.y,atom2.z}}
			table.insert(bonds,bond)
			end
		end
	end
end
	if res["C  "]~=nil and i<#box then if box[i+1]["N  "]~=nil then
	if math.sqrt((res["C  "].x-box[i+1]["N  "].x)^2+(res["C  "].y-box[i+1]["N  "].y)^2+(res["C  "].z-box[i+1]["N  "].z)^2) < maxbondlen then
	peptide={{res["C  "].x,res["C  "].y,res["C  "].z},{box[i+1]["N  "].x,box[i+1]["N  "].y,box[i+1]["N  "].z}}
	table.insert(bonds,peptide)
	end end
	end
	if res["O3'"]~=nil and i<#box then if box[i+1]["P  "]~=nil then
	if math.sqrt((res["O3'"].x-box[i+1]["P  "].x)^2+(res["O3'"].y-box[i+1]["P  "].y)^2+(res["O3'"].z-box[i+1]["P  "].z)^2) < maxbondlen then
	phospho={{res["O3'"].x,res["O3'"].y,res["O3'"].z},{box[i+1]["P  "].x,box[i+1]["P  "].y,box[i+1]["P  "].z}}
	table.insert(bonds,phospho)
	end end
	end
if #bonds>0 then res["bonds"]=bonds	end
end  --##


if thing==protein then drawbonds(box) else drawbonds(box,true) end -- draws the normal aa bonds, true for symm_box


nextra=#extrabox   -- gets SS bonds
if nextra>1 then
    extrabonds={}
    for j,atom in ipairs(extrabox) do
	    for m,n in pairs(atom) do at1=m end  --get atom key
	    for k,atom2 in ipairs(extrabox) do
		    if k>j then
			    for o,p in pairs(atom2) do at2=o end  --ditto
			    d=math.sqrt((extrabox[j][at1].x-extrabox[k][at2].x)^2+(extrabox[j][at1].y-extrabox[k][at2].y)^2+(extrabox[j][at1].z-extrabox[k][at2].z)^2)
			    if d<maxssbondlen then
			    bond={{extrabox[j][at1].x,extrabox[j][at1].y,extrabox[j][at1].z},{extrabox[k][at2].x,extrabox[k][at2].y,extrabox[k][at2].z}}
			    table.insert(extrabonds,bond)
			    end
		    end
	    end
    end
    extraboxbond={{["ssbrs"]={},["bonds"]=extrabonds}}  --dummy residue containing SS bonds
    drawbonds(extraboxbond,true)  -- draws any SSBRs
end
--for i,j in ipairs(extraboxbond) do for k,l in pairs(j) do print(k) if k=="bonds" then for m,n in ipairs(l) do print(n[1][1],n[1][2],n[1][3],n[2][1],n[2][2],n[2][3]) end end end end
end


function getBoxAtoms()  --££
group = js.new(threejs.Group)
if rotamat~=nil then group:applyMatrix4(rotamat) end
getCentre()
getBonds(protein)
drawatoms()
if gen_symm==true then
   if symm_box_filled~=true then get_centroid() end
   getBonds(symm_box)
   drawatoms()
end
end  --££


function drawbonds(abox,alt)
if alt==true then matter=altmat alt=false else matter=material end
vector1=js.new(threejs.Vector3,0,0,0)
vector2=js.new(threejs.Vector3,0,0,0)
geobond = js.new(threejs.BufferGeometry)	
--for i,res in ipairs(protein) do
for i,res in ipairs(abox) do
	if res["bonds"]~=nil then
	for j,bond in ipairs(res["bonds"]) do 
        vector1:set(bond[1][1]-xcent,bond[1][2]-ycent,bond[1][3]-zcent)
        vector2:set(bond[2][1]-xcent,bond[2][2]-ycent,bond[2][3]-zcent)
		vectors=window:Array(vector1,vector2)
		geometry1=geobond:clone()
		geometry1:setFromPoints(vectors)
		mesh1 = js.new(threejs.Line, geometry1, matter)
		group:add(mesh1)
	end
	end
end
--for i,res in ipairs(box) do print("Residue:",i) for j,k in pairs(res) do print(j) end end
end


function drawatoms() 
vectoratom=js.new(threejs.Vector3,0,0,0)
atgeom=js.new(threejs.BufferGeometry)
for i,j in ipairs(box) do 
		for k,l in pairs(j) do
			attype=k:sub(1,1) 
			--if attype=="N" or attype=="O" or attype=="S" or attype=="P" then
                        if attype~="C" and attype~="b" and attype~="s" then  --"b" and "s" are the dummy atom types for the bonds and SSBR's
                        if attype=="N" then atcol="blue" elseif attype=="O" then atcol="red" elseif attype=="S" then atcol="#0FFF50" else atcol="#0FFF50" end
				--print(attype,atcol,l.x-xcent,l.y-ycent,l.z-zcent)
				materat=js.new(threejs.PointsMaterial)
				materat.size=3
				materat.sizeAttenuation=false
				materat.color:set(atcol)
				vectorat=window:Array(vectoratom:set(l.x-xcent,l.y-ycent,l.z-zcent))
			geomat=atgeom:clone()
			geomat:setFromPoints(vectorat)
			meshat=js.new(threejs.Points,geomat,materat)
			group:add(meshat)
			end
		end
end
--Next lines not needed in map version
scene:add(group)
if spin==true then renderer:setAnimationLoop(twirl)	else renderer:setAnimationLoop() renderer:render(scene, camera) end
document:getElementById("mutate").disabled=false
document:getElementById("move").disabled=false
end


function twirl()
	theta=0.02
        twang=twang+theta
        group:rotateOnWorldAxis(xworld,theta*math.cos(twang))
        group:rotateOnWorldAxis(yworld,theta*math.sin(twang))
	renderer:render(scene, camera)
end


function uploadMTZ()
document:getElementById("mtztips").innerHTML="Please wait for uploads to complete!" 
--document:getElementById("mapdownload").disabled=true
diffneeded=document:getElementById("diffmap").checked
formData=js.new(window.FormData)
formData:append("pdbfile",pdbfile)
formData:append("mtzfile",mapfile)
formData:append("diffneeded",diffneeded)
fetchObject=js.new(window.FormData)
fetchObject.method="POST"
fetchObject.body=formData
fetchpromise=window:fetch("/cgi-bin/run_gemmi.py",fetchObject)
fetchpromise["then"](fetchpromise, getMapBack)  --override lua's own if/then command
end


function getMapBack()
document:getElementById("mtztips").innerHTML="Please wait for MAP file to download."
if diffneeded~=true then getmap=nameparts[1]..".map" else getmap=nameparts[1].."_diff.map" end
mapURI="/temp/"..getmap
mapdownload=document:createElement("a")
mapdownload.href=mapURI
mapdownload.download=getmap
mapdownload:click()
document:getElementById("mtztips").innerHTML="Please use map file chooser button to load it."
end


function readmapfile()
mapfile=mapinput.files[0]
name1=pdbfile.name name2=mapfile.name print(name1, name2)
nameparts={} for word in string.gmatch(mapfile.name,"[^.]+") do table.insert(nameparts,word) end --print(nameparts[1],nameparts[2])
ref_filetype=string.lower(nameparts[2])
if ref_filetype=="mtz" and nAt>0 then
	--mtztext="Tick for difference map: <input type=\"checkbox\" id=\"diffmap\">&nbsp; &nbsp;<button id=\"uploadmtz\">Make map from MTZ file</button>&nbsp; &nbsp;<button id=\"mapdownload\" disabled>Download map file</button> <p> <small><span id=\"mtztips\"></span></small>"
        mtztext="Tick for difference map: <input type=\"checkbox\" id=\"diffmap\">&nbsp; &nbsp;<button id=\"uploadmtz\">Make map from MTZ file</button><p><small><span id=\"mtztips\"></span></small>"
        document:getElementById("info").innerHTML=mtztext
	mtzin=document:getElementById("uploadmtz")
	mtzin:addEventListener("click", uploadMTZ)
else	
	mapreader=js.new(window.FileReader)
	mapreader.onload=readmap
	mapreader:readAsArrayBuffer(mapfile)
	document:getElementById("info").innerHTML="Please wait"
end
end


function readmap()
map_fills_cell=false
map=mapreader.result
viewer=js.new(window.DataView,map)
print("Grid:")
grid={} for g=1,10 do grid[g]=viewer:getInt32((g-1)*4,true) print(grid[g]) end	--up to 40. Grid is num col, row, sect, mode, start col, row, sect, num div on a, b, c. Mode is 2 for 32 bit floats
na=grid[8] nb=grid[9] nc=grid[10]
print("Cell:")
cell={} for c=1,6 do cell[c]=viewer:getFloat32(40+(c-1)*4,true) print(cell[c]) end	--up to 64
print("Axis order:")
order={} for o=1,3 do order[o]=viewer:getInt32(64+(o-1)*4,true) print(order[o]) end	--2,1,3 z sections (ncode 1), 3,1,2 y sections (ncode3), 1,2,3 for gemmi, 3,2,1 for EM
print("Min, max, mean density:")
mimameden={} for m=1,3 do mimameden[m]=viewer:getFloat32(76+(m-1)*4,true) print(mimameden[m]) end mean=mimameden[3]
print("Symm:")
sym={} for s=1,3 do sym[s]=viewer:getInt32(88+(s-1)*4,true) print(sym[s]) end denoffset=1024+sym[2] -- sym[1]=sg num, sym[2] is length of string of symops (not used)
rms=viewer:getFloat32(216,true) print("rms:",rms)
nlbl=viewer:getInt32(220,true) --print(nlbl)

lblbufarr=map:slice(224,1023)
lblvwr=js.new(window.TextDecoder,"utf-8")
labels=lblvwr:decode(lblbufarr) --print(labels) 

symbufarr=map:slice(1024,denoffset);
symvwr=js.new(window.TextDecoder,"utf-8")
symops=symvwr:decode(symbufarr) --print(symops) 	
	
--den={} for d=1,6 do den[d]=viewer:getFloat32(denoffset+(d-1)*4,true) print(den[d]) end

if order[1]==2 and order[2]==1 and order[3]==3 then print("Old normal CCP4 map") ncode=1
amin=grid[6] numona=grid[2] bmin=grid[5] numonb=grid[1] cmin=grid[7] numonc=grid[3]
elseif order[1]==3 and order[2]==1 and order[3]==2 then print("Old Mono- or tri-clinic CCP4 map") ncode=3
amin=grid[6] numona=grid[2] bmin=grid[7] numonb=grid[3] cmin=grid[5] numonc=grid[1]
elseif order[1]==1 and order[2]==2 and order[3]==3 then print("Gemmi map") ncode="gemmi"
amin=grid[5] numona=grid[1] bmin=grid[6] numonb=grid[2] cmin=grid[7] numonc=grid[3]
elseif order[1]==3 and order[2]==2 and order[3]==1 then print("CCPEM map") ncode="ccpem"
amin=grid[7] numona=grid[3] bmin=grid[6] numonb=grid[2] cmin=grid[5] numonc=grid[1]
else document:getElementById("info").innerHTML="Bad fast medium slow axis order" return
end
a,b,c,alpha_deg,beta_deg,gamma_deg=cell[1],cell[2],cell[3],cell[4],cell[5],cell[6]

if amin==0 and bmin==0 and cmin==0 and numona==na and numonb==nb and numonc==nc then map_fills_cell=true print("Map covers unit cell rather than molecule") end  --occurs with some EBI PDBe maps

fractorth()
getMapBox()
end	


function getMapBox()	
document:getElementById("info").innerHTML=centrat.res.." "..centrat.num..centrat.chain
if foundSideCen==false then
   if protein[startpos]["CB "]~=nil then mcenat=protein[startpos]["CB "] elseif protein[startpos]["CA "]~=nil then mcenat=protein[startpos]["CA "] else mcenat=centrat end
else
   mcenat=centrat
   foundSideCen=false
end
print("Map centre:",mcenat.x,mcenat.y,mcenat.z)
xmin=mcenat.x-maprad xmax=mcenat.x+maprad ymin=mcenat.y-maprad ymax=mcenat.y+maprad zmin=mcenat.z-maprad zmax=mcenat.z+maprad
	
clo=window:Array(xmin,ymin,zmin)	
far=window:Array(xmax,ymax,zmax)
fclo=mathjs:multiply(Linv,clo)
ffar=mathjs:multiply(Linv,far)
--print(fclo,ffar)	

xminf=(mathjs:subset(fclo,el1))	yminf=(mathjs:subset(fclo,el2)) zminf=(mathjs:subset(fclo,el3))	
xmaxf=(mathjs:subset(ffar,el1))	ymaxf=(mathjs:subset(ffar,el2)) zmaxf=(mathjs:subset(ffar,el3))	
--print(xminf,xmaxf,yminf,ymaxf,zminf,zmaxf)
	
mapamin=math.ceil(xminf*na) intamin=mapamin-amin mapamax=math.ceil(xmaxf*na) intamax=mapamax-amin ndiva=(mapamax-mapamin) --math.abs maybe needed!! for -ve on axis
mapbmin=math.ceil(yminf*nb) intbmin=mapbmin-bmin mapbmax=math.ceil(ymaxf*nb) intbmax=mapbmax-bmin ndivb=(mapbmax-mapbmin)
mapcmin=math.ceil(zminf*nc) intcmin=mapcmin-cmin mapcmax=math.ceil(zmaxf*nc) intcmax=mapcmax-cmin ndivc=(mapcmax-mapcmin)
--print(intamin,intamax,intbmin,intbmax,intcmin,intcmax,"ndivs",ndiva,ndivb,ndivc)

minimap={}
for  cc=intcmin,intcmax do
			if map_fills_cell==true then cc=cc%nc if cc<0 then cc=nc-cc end end
			minimapsection={}
			for aa=intamin,intamax do
				if map_fills_cell==true then aa=aa%na if aa<0 then aa=na-aa end end
				minimaprow={}
				for bb=intbmin,intbmax do
					if map_fills_cell==true then bb=bb%nb if bb<0 then bb=nb-bb end end
					if ncode==1 then
						posi=denoffset+cc*numona*numonb*4+aa*numonb*4+bb*4
					elseif ncode==3 then
						posi=denoffset+bb*numona*numonc*4+aa*numonc*4+cc*4
					elseif ncode=="gemmi" then
						posi=denoffset+cc*numona*numonb*4+bb*numona*4+aa*4
					elseif ncode=="ccpem" then
						posi=denoffset+aa*numonc*numonb*4+bb*numonc*4+cc*4
					end
					dens=viewer:getFloat32(posi,true)
					normdens=(dens-mean)/rms
					table.insert(minimaprow,normdens)
end table.insert(minimapsection,minimaprow) end table.insert(minimap,minimapsection) end

--for i,j in ipairs(minimap) do for k,l in ipairs(j) do for m,n in ipairs(l) do print(n)	end end end 
	
ndiv=window:Array(ndiva,ndivb,ndivc)	
start=window:Array(0,0,0)
bounds=window:Array(start,ndiv)
isomesh=window:surfaceNets(ndiv, function(_,x,y,z) return denread(x,y,z) end, bounds)
window.console:log(isomesh)
posis=isomesh.positions	

vertices={}	
	
for i=0,#posis-1 do
	xcontf=(mapamin+posis[i][0])/na  ycontf=(mapbmin+posis[i][1])/nb  zcontf=(mapcmin+posis[i][2])/nc 
	minimapxyzf=window:Array(xcontf,ycontf,zcontf)
	minimapxyzo=mathjs:multiply(L,minimapxyzf)
	minixo=mathjs:subset(minimapxyzo,el1) miniyo=mathjs:subset(minimapxyzo,el2) minizo=mathjs:subset(minimapxyzo,el3)
	table.insert(vertices,{minixo-xcent,miniyo-ycent,minizo-zcent})
end

--scene:add(group)

vector1=js.new(threejs.Vector3,0,0,0)
vector2=js.new(threejs.Vector3,0,0,0)
vector3=js.new(threejs.Vector3,0,0,0)
geomap = js.new(threejs.BufferGeometry)

triangs=isomesh.cells		
for i=0,#triangs-1 do
	verts=triangs[i]
	vert1=verts[0]+1 vert2=verts[1]+1 vert3=verts[2]+1  -- +1 to convert to lua indices in vertices table	
	--print("\n",vert1,vert2,vert3)
	vector1:set(vertices[vert1][1],vertices[vert1][2],vertices[vert1][3])
	vector2:set(vertices[vert2][1],vertices[vert2][2],vertices[vert2][3])
	vector3:set(vertices[vert3][1],vertices[vert3][2],vertices[vert3][3])
	vectors=window:Array(vector1,vector2,vector3)
	mapgeo=geomap:clone()
	mapgeo:setFromPoints(vectors)
	meshmap = js.new(threejs.LineLoop, mapgeo, matermap)
	group:add(meshmap)	
end	
scene:add(group)
--if spin==true then renderer:setAnimationLoop(function() twirl() end)	else offsetX,offsetY=0,0 renderer:setAnimationLoop() renderer:render(scene, camera) document:getElementById("move").disabled=false end
if spin==true then renderer:setAnimationLoop(twirl)	else renderer:setAnimationLoop() renderer:render(scene, camera) end
document:getElementById("mutate").disabled=false
document:getElementById("move").disabled=false

geobond:dispose() geometry1:dispose() atgeom:dispose() geomat:dispose() geomap:dispose() --RAM housekeeping
material:dispose() altmat:dispose() matermap:dispose() matter:dispose() materat:dispose()

end


--Dictionaries
dictable={   --target distances are read from dist table and are grouped to be within about 0.1 Angstroms. *prev/next*=ignored restraints
["N  "]={["CA "]=4,["CB "]=9,["C  "]=9}, --*prev*["C  "],*prev*["CA "],*prev*["O  "]
["CA "]={["N  "]=4,["C  "]=4,["CB "]=4,["O  "]=8,["CG "]=9,["OG "]=8,["OG1"]=8,["CG1"]=9,["CG2"]=9,["SG "]=10}, --*prev*["C  "],*next*["N  "]
["C  "]={["CA "]=4,["CB "]=9,["N  "]=9,["O  "]=1}, --*next*["N  "],*next*["CA "],*next pro*["CD "]
["O  "]={["C  "]=1,["CA "]=8}, --*next*["N  "]
["CB "]={["CA "]=4,["N  "]=9,["C  "]=9,["CG "]=4,["OG "]=3,["OG1"]=3,["CG1"]=4,["CG2"]=4,["SG "]=5,["CD "]=9,["CD1"]=9,["CD2"]=9,["OD1"]=8,["OD2"]=8,["ND2"]=8,["ND1"]=9,["SD "]=10},
["CG "]={["CA "]=9,["CB "]=4,["CD1"]=4,["CD2"]=4,["CD "]=4,["SD "]=5,["OD1"]=1,["OD2"]=1,["ND2"]=2,["ND1"]=2,["OE1"]=8,["OE2"]=8,["NE1"]=6,["NE2"]=8,["NE "]=9,["CE1"]=8,["CE2"]=8,["CE3"]=9,["CE "]=10}, --Met/Lys CG/CE sorted in code
["OG "]={["CB "]=3,["CA "]=8},
["SG "]={["CB "]=5,["CA "]=10},
["OG1"]={["CB "]=3,["CA "]=8,["CG2"]=8},
["CG1"]={["CB "]=4,["CD1"]=4,["CA "]=9,["CG2"]=9},
["CG2"]={["CB "]=4,["CA "]=9,["OG1"]=8,["CG1"]=9},
["CD "]={["CG "]=4,["CB "]=9,["OE1"]=1,["OE2"]=1,["NE2"]=2,["CE "]=4,["NE "]=4,["CZ "]=9,["NZ "]=9}, --*pro prev*["C  "]
["CD1"]={["CG "]=4,["CG1"]=4,["CB "]=9,["CZ "]=8,["CD2"]=8,["CE1"]=3,["NE1"]=3},   --cd1/2 aliphatic/aromatic simplified
["CD2"]={["CG "]=4,["CE2"]=3,["CB "]=9,["CZ "]=8,["CE3"]=3,["NE1"]=6,["CZ2"]=9,["CZ3"]=8,["NE2"]=3,["ND1"]=6,["CD1"]=8},
["OD1"]={["CG "]=1,["CB "]=8,["OD2"]=6,["ND2"]=6},
["OD2"]={["CG "]=1,["CB "]=8,["OD1"]=6},
["ND1"]={["CG "]=3,["CB "]=9,["CE1"]=2,["CD2"]=6,["NE2"]=6},
["ND2"]={["CG "]=2,["CB "]=8,["OD1"]=6},
["SD "]={["CG "]=5,["CE "]=5,["CB "]=10},
["CE "]={["SD "]=5,["NZ "]=4,["CD "]=4,["CG "]=10},
["CE1"]={["CD1"]=3,["CG "]=8,["CZ "]=3,["OH "]=8,["ND1"]=2,["NE2"]=2,["CE2"]=8},  --ce1/2 same as cd1/2
["CE2"]={["CD2"]=3,["CG "]=8,["CZ "]=3,["OH "]=8,["NE1"]=3,["CZ2"]=3,["CH2"]=8,["CE3"]=8,["CE1"]=8},
["CE3"]={["CD2"]=3,["CG "]=9,["CZ3"]=3,["CH2"]=8,["CE2"]=8},
["OE1"]={["CD "]=1,["CG "]=8,["OE2"]=6,["NE2"]=6},
["OE2"]={["CD "]=1,["CG "]=8,["OE1"]=6},
["NE "]={["CD "]=4,["CG "]=9,["CZ "]=2,["NH1"]=7,["NH2"]=7},
["NE1"]={["CD1"]=3,["CG "]=6,["CE2"]=3,["CD2"]=6,["CZ2"]=9},
["NE2"]={["CD "]=2,["CG "]=8,["OE1"]=6,["CD2"]=2,["CE1"]=2,["ND1"]=6},
["CZ "]={["CE1"]=3,["CE2"]=3,["CD1"]=8,["CD2"]=8,["OH "]=3,["NE "]=2,["CD "]=9,["NH1"]=2,["NH2"]=2},
["NZ "]={["CE "]=4,["CD "]=9},
["CZ2"]={["CE2"]=3,["CD2"]=9,["NE1"]=9,["CH2"]=3,["CZ3"]=8},
["CZ3"]={["CE3"]=3,["CH2"]=3,["CD2"]=8,["CZ2"]=8},
["CH2"]={["CZ2"]=3,["CZ3"]=3,["CE2"]=8,["CE3"]=8},
["NH1"]={["CZ "]=2,["NH2"]=7,["NE "]=7},
["NH2"]={["CZ "]=2,["NH1"]=7,["NE "]=7},
["OH "]={["CZ "]=3,["CE1"]=8,["CE2"]=8}
}

protable={   --target distances for proline
["N  "]={["CA "]=4,["CB "]=8,["C  "]=9,["CG "]=8,["CD "]=4}, --*prev*["C  "],*prev*["CA "],*prev*["O  "]
["CA "]={["N  "]=4,["C  "]=4,["CB "]=4,["O  "]=8,["CG "]=8,["CD "]=9}, --*prev*["C  "],*next*["N  "]
["C  "]={["CA "]=4,["CB "]=9,["N  "]=9,["O  "]=1}, --*next*["N  "],*next*["CA "],*next pro*["CD "]
["O  "]={["C  "]=1,["CA "]=8}, --*next*["N  "]
["CB "]={["CA "]=4,["N  "]=8,["C  "]=9,["CG "]=4,["CD "]=9},
["CG "]={["CA "]=8,["CB "]=4,["CD "]=4,["N  "]=8},
["CD "]={["CG "]=4,["CB "]=9,["N  "]=4,["CA "]=9} --*pro prev*["C  "]
}

tetrigtable={  --tetragonal or trigonal carbon centres, neighbours and distance from centroid of neighbours. CA and Thr/Ile CB chirality not checked.
["ASP"]={["CG "]={"CB ","OD1","OD2",0.1}},["ASN"]={["CG "]={"CB ","OD1","ND2",0.1}},
["GLU"]={["CD "]={"CG ","OE1","OE2",0.1}},["GLN"]={["CD "]={"CG ","OE1","NE2",0.1}},
["PHE"]={["CG "]={"CB ","CD1","CD2",0.0}},
["TYR"]={["CG "]={"CB ","CD1","CD2",0.0},["CZ "]={"CE1","CE2","OH ",0.0}},
["TRP"]={["CG "]={"CB ","CD1","CD2",0.07},["CD2"]={"CG ","CE2","CE3",0.1},["CE2"]={"CD2","NE1","CZ2",0.1}},
["LEU"]={["CG "]={"CB ","CD1","CD2",0.5}},["VAL"]={["CB "]={"CA ","CG1","CG2",0.5}},
["ILE"]={["CB "]={"CA ","CG1","CG2",0.5}},["THR"]={["CB "]={"CA ","OG1","CG2",0.5}},
["ARG"]={["CZ "]={"NE ","NH1","NH2",0.0}},["HIS"]={["CG "]={"CB ","ND1","CD2",0.04}}}

planetable={  --restrain some interatomic vectors to be parallel e.g. CD1-CE1 bond is parallel to CB-CG in Phe/Tyr but shorter by factor of 0.91.
["PHE"]={["CE1"]={"CB ","CG ","CD1",0.91},["CE2"]={"CB ","CG ","CD2",0.91},["CZ "]={"CG ","CD1","CE2",1.0},["CZ "]={"CG ","CD2","CE1",1.0}},
["TYR"]={["CE1"]={"CB ","CG ","CD1",0.91},["CE2"]={"CB ","CG ","CD2",0.91},["CZ "]={"CG ","CD1","CE2",1.0},["CZ "]={"CG ","CD2","CE1",1.0},["OH "]={"CB ","CG ","CZ ",0.93}},
["TRP"]={["CH2"]={"CD2","CE2","CZ3",1.0},["CH2"]={"CD2","CE3","CZ2",1.0},["CZ3"]={"CE2","CZ2","CE3",1.0},["CZ2"]={"CG ","NE1","CE3",1.26},["CH2"]={"CG ","NE1","CZ3",0.63}},
["ARG"]={["NH2"]={"CD ","NE ","CZ ",1.0}}}


function addTetrig()  --calculate distance of central atom from centroid of neighbours and the diff from ideal
ctat=residue[centat] nbr1=residue[nayb1]  nbr2=residue[nayb2]  nbr3=residue[nayb3]
xcentet=(nbr1.x+nbr2.x+nbr3.x)/3
ycentet=(nbr1.y+nbr2.y+nbr3.y)/3
zcentet=(nbr1.z+nbr2.z+nbr3.z)/3
dx=ctat.x-xcentet dy=ctat.y-ycentet dz=ctat.z-zcentet
dcalc=math.sqrt(dx^2+dy^2+dz^2)
diff=dtarg-dcalc
if dcalc~=0 then
table.insert(atDeltas,{["diff"]=diff,["at1"]="tetrig",["at2"]=centat,["delx"]=diff*dx/dcalc,["dely"]=diff*dy/dcalc,["delz"]=diff*dz/dcalc})
end
--print("In addTetrig for:", centat, nayb1, nayb2, nayb3,"dcalc:",dcalc)
--print(ctat.x,ctat.y,ctat.z,nbr1.x,nbr1.y,nbr1.z,nbr2.x,nbr2.y,nbr2.z,nbr3.x,nbr3.y,nbr3.z)
end

function addPlane()  --calculate distance of atom from ideal position that would make the vectors parallel
planx=residue[tail2].x+(residue[tip1].x-residue[tail1].x)*dscale
plany=residue[tail2].y+(residue[tip1].y-residue[tail1].y)*dscale
planz=residue[tail2].z+(residue[tip1].z-residue[tail1].z)*dscale
dx=residue[planat].x-planx dy=residue[planat].y-plany dz=residue[planat].z-planz
dcalc=math.sqrt(dx^2+dy^2+dz^2)
diff=dtarg-dcalc
if dcalc~=0 then
table.insert(atDeltas,{["diff"]=diff,["at1"]="plane",["at2"]=planat,["delx"]=diff*dx/dcalc,["dely"]=diff*dy/dcalc,["delz"]=diff*dz/dcalc})
end
--print("In addPlane")
end

function addTrpPlane()  --calculate distance of CD1 from its ideal position on the line bisecting CH2-CZ3 and CD2-CE2 bonds. Helps to make the Trp flat.
dum1x=(residue["CZ3"].x+residue["CH2"].x)/2 dum1y=(residue["CZ3"].y+residue["CH2"].y)/2 dum1z=(residue["CZ3"].z+residue["CH2"].z)/2
dum2x=(residue["CD2"].x+residue["CE2"].x)/2 dum2y=(residue["CD2"].y+residue["CE2"].y)/2 dum2z=(residue["CD2"].z+residue["CE2"].z)/2
planx=dum1x+(dum2x-dum1x)*1.89 plany=dum1y+(dum2y-dum1y)*1.89 planz=dum1z+(dum2z-dum1z)*1.89
dx=residue["CD1"].x-planx dy=residue["CD1"].y-plany dz=residue["CD1"].z-planz
dcalc=math.sqrt(dx^2+dy^2+dz^2)
diff=dtarg-dcalc
if dcalc~=0 then
table.insert(atDeltas,{["diff"]=diff,["at1"]="plane",["at2"]="CD1",["delx"]=diff*dx/dcalc,["dely"]=diff*dy/dcalc,["delz"]=diff*dz/dcalc})
end
--print("In addTrpPlane")
end

function checkChir()  --calculate the chiral volume
chrat=residue[chirat] nbr1=residue[nayb1]  nbr2=residue[nayb2]  nbr3=residue[nayb3]
detRow1=window:Array(nbr1.x-chrat.x, nbr1.y-chrat.y, nbr1.z-chrat.z)
detRow2=window:Array(nbr2.x-chrat.x, nbr2.y-chrat.y, nbr2.z-chrat.z)
detRow3=window:Array(nbr3.x-chrat.x, nbr3.y-chrat.y, nbr3.z-chrat.z)
chirDet=window:Array(detRow1,detRow2,detRow3)
chivdet=js.new(mathjs.matrix,chirDet)
chiv=mathjs:det(chivdet)
print("In checkChir Chiv:",chiv,"for:",chirat,nayb1,nayb2,nayb3)
--print(chrat.x,chrat.y,chrat.z,nbr1.x,nbr1.y,nbr1.z,nbr2.x,nbr2.y,nbr2.z,nbr3.x,nbr3.y,nbr3.z)
--if chiv<0 and math.abs(chiv)>0.1 then correctChir() end    --invert only if |chiv| is big enough to avoid divide by zero error
if chiv<0 then correctChir() end
end

function correctChir()
xcentet=(nbr1.x+nbr2.x+nbr3.x)/3  ycentet=(nbr1.y+nbr2.y+nbr3.y)/3  zcentet=(nbr1.z+nbr2.z+nbr3.z)/3
mov1=chrat  mov2=nbr3
dx=mov1.x-xcentet dy=mov1.y-ycentet dz=mov1.z-zcentet  dcalc=math.sqrt(dx^2+dy^2+dz^2)  diff=-0.5 --diff=dtarg-dcalc
if dcalc~=0 then
delx=diff*dx/dcalc  dely=diff*dy/dcalc  delz=diff*dz/dcalc
mov1.x=mov1.x+delx  mov1.y=mov1.y+dely  mov1.z=mov1.z+delz
mov2.x=mov2.x-delx  mov2.y=mov2.y-dely  mov2.z=mov2.z-delz
end
--print("In correctChir")
end


function mmm()    --massively minimal minimiser (mmm) mainly for side chains (pat not pending)  

if ncyc>0 and ncyc%20==0 and resnam~="GLY" then chirat="CA " nayb1="N  " nayb2="C  " nayb3="CB " checkChir() end  --Check the main chain chiral volume at intervals
if resnam=="THR" or resnam=="ILE" then  --Check the side chain chiral volume at bigger intervals
  if ncyc>0 and ncyc%33==0 then
    chirAtObjs=tetrigtable[resnam]
    for k,l in pairs(chirAtObjs) do  --$$
        --print(k,l)
        chirat=k nayb1=l[1] nayb2=l[2] nayb3=l[3] checkChir()
    end  --$$
  end
end

dist={1.2,1.3,1.4,1.5,1.8,2.2,2.3,2.4,2.5,2.8}
--     1   2   3   4   5   6   7   8   9  10   simplified distances
atDeltas={}
--for k,l in ipairs(lookup1) do print(l) end
kk=0
for k,l in pairs(residue) do
kk=kk+1  mm=0
for m,n in pairs(residue) do
mm=mm+1

if mm>kk then   --calculate the deviations from target distances (diff)
      --print(k,l.x,l.y,l.z)
      if dictable[k]~=nil and dictable[m]~=nil then
          if resnam~="PRO" then dtarg=dist[dictable[k][m]] else dtarg=dist[protable[k][m]] end
          if resnam=="LYS" then
					    if k=="CG " and m=="CE " then dtarg=2.5
					    elseif k=="CE " and m=="CG " then dtarg=2.5
					    end
          end
          if dtarg~=nil then  --print(k,m)
                 dx=n.x-l.x
                 dy=n.y-l.y
                 dz=n.z-l.z
                 dcalc=math.sqrt(dx^2+dy^2+dz^2)
                 diff=dtarg-dcalc
                 table.insert(atDeltas,{["diff"]=diff,["at1"]=k,["at2"]=m})
          end
      end
end

end
end

--Calculate the tetragonal and trigonal target distances
if resnam~="GLY" then centat="CA " nayb1="N  " nayb2="CB " nayb3="C  " dtarg=0.5 addTetrig() end
if tetrigtable[resnam]~=nil then  --@@
    centAtObjs=tetrigtable[resnam]
    for k,l in pairs(centAtObjs) do  --$$
        centat=k nayb1=l[1] nayb2=l[2] nayb3=l[3] dtarg=l[4] addTetrig()
    end  --$$
end  --@@

if ncyc>0 and maxdiff<0.3 and ncyc%3==0 then  --trickle in the planarity targets
if planetable[resnam]~=nil then  --@@
    planeAtObjs=planetable[resnam]
    for k,l in pairs(planeAtObjs) do  --$$
        planat=k tail1=l[1] tip1=l[2] tail2=l[3] dscale=l[4] addPlane()
    end  --$$
end  --@@
if resnam=="TRP" then addTrpPlane() end --extra planarity target for CD1 of Trp
end

ncyc=ncyc+1
--print("Cycle number:",ncyc)
--for i,j in ipairs(atDeltas) do for k,l in pairs(j) do print(l) end end
maxdiff=-100
for i,j in ipairs(atDeltas) do absdiff=math.abs(j.diff)  --find the worst distance violation
    if absdiff>maxdiff then
        maxdiff=absdiff
        worstPair={j.at1,j.at2}
        if j.at1=="tetrig" or j.at1=="plane" then delx=j.delx dely=j.dely delz=j.delz
        end
    end
end

if ncyc>=100 and maxdiff<0.1 then sorted=true end
--for i,j in ipairs(worstPair) do print(i,j) end
at1=worstPair[1]
at2=worstPair[2]
--print("The worst pair (atom name) with delta:",at1,at2,maxdiff)
if at1=="tetrig" or at1=="plane" then --##
mov=residue[at2]
--print(at1)
else
--print("d")
sporad1=math.random(1,2)  --decide at random which atom to move for atom-atom distance violations.
if sporad1==1 then fixed=at1 moving=at2 else fixed=at2 moving=at1 end
if moving=="N  " then moving=fixed fixed="N  " end  --main chain N and C fixed
if moving=="C  " then moving=fixed fixed="C  " end
--print("final decision fixed:",fixed,"moving:",moving)
fix=residue[fixed]
mov=residue[moving]
dx=mov.x-fix.x
dy=mov.y-fix.y
dz=mov.z-fix.z
dcalc=math.sqrt(dx^2+dy^2+dz^2)
if resnam~="PRO" then dtarg=dist[dictable[fixed][moving]] else dtarg=dist[protable[fixed][moving]] end
if resnam=="LYS" then 
	if fixed=="CG " and moving=="CE " then dtarg=2.5
	elseif fixed=="CE " and moving=="CG " then dtarg=2.5
	end
end		
diff=dtarg-dcalc
delx=diff*dx/dcalc dely=diff*dy/dcalc delz=diff*dz/dcalc
end  --##

mov.x=mov.x+delx/2  mov.y=mov.y+dely/2  mov.z=mov.z+delz/2  -- Dampen the shifts
--mov.x=mov.x+delx  mov.y=mov.y+dely  mov.z=mov.z+delz

end  --mmm


function tidy()
document:getElementById("info").innerHTML="Please wait"
window:setTimeout(function() tidy2() end, 10)  --10 msec delay to let the DOM display the message
end

function tidy2()
residue=protein[startpos]
resnam=string.sub(residue["CA "].res,1,3)
ncyc=0   --normally needed
sorted=false
--mmm() --normally commented out
--repeat mmm() until sorted or ncyc>400
repeat mmm() until sorted or ncyc>200
print("Residue:",resnam,startpos,"max delta:",maxdiff,"after", ncyc, "cycles.")
rerun() --needed normally
end


function denread(xiso,yiso,ziso)
--print("denread:",xiso,yiso,ziso)	
return minimap[ziso+1][xiso+1][yiso+1]-contour_level		
end	


function getStartPos()  --if startpos < 1st resnum go to 1st res in chain, if startpos > last resnum go to last res in chain
for i,res in ipairs(protein) do
    _,firstAtom=next(res) --print(firstAtom.chain,firstAtom.num)
    firstAtomResNum=string.sub(firstAtom.num,1,-2)  firstAtomResNum=tonumber(firstAtomResNum) --remove insertion code so resnum is a number
    if chain==firstAtom.chain then
       lastInChain=i
       if startpos<=firstAtomResNum then wantedpos=i break else wantedpos=lastInChain end
    end
end
startpos=wantedpos  print("Actual res num, index:",firstAtomResNum,wantedpos)
position.value=firstAtomResNum  --send resnum to the form ignoring insertion code
chainbox.value=firstAtom.chain  --send chain to the form
end


function readform()
document:getElementById("info").innerHTML="Please wait"
startpos=tonumber(position.value)
chainbox=document:getElementById("chain")
if chainbox~=js.null then chain=chainbox.value chain=string.upper(chain) chain=string.format("%2s",chain) getStartPos() end
contour=document:getElementById("contour")
contour_level=tonumber(contour.value) if contour_level<0 then mapcol=mapcolneg else mapcol=mapcolpos end 
matermap.color:set(mapcol)
moving=false
if document:getElementById("symmetry").checked==true then gen_symm=true else gen_symm=false end
end


function swapdensign()
contour_level=tonumber(document:getElementById("contour").value)
contour_level=-contour_level
document:getElementById("contour").value=contour_level
--print(document:getElementById("contour").value)
if minimap~=nil then if #(minimap)>0 then run() end end	
end


sidechain_atoms={
["ALA"]={"CB "},
["ARG"]={"CB ", "CG ", "CD ", "NE ", "CZ ", "NH1", "NH2"},
["ASN"]={"CB ", "CG ", "OD1", "ND2"},
["ASP"]={"CB ", "CG ", "OD1", "OD2"},
["CYS"]={"CB ", "SG "},
["GLN"]={"CB ", "CG ", "CD ", "OE1", "NE2"},
["GLU"]={"CB ", "CG ", "CD ", "OE1", "OE2"},
["GLY"]={},
["HIS"]={"CB ", "CG ", "ND1", "CD2", "CE1", "NE2"},
["ILE"]={"CB ", "CG1", "CG2", "CD1"},
["LEU"]={"CB ", "CG ", "CD1", "CD2"},
["LYS"]={"CB ", "CG ", "CD ", "CE ", "NZ "},
["MET"]={"CB ", "CG ", "SD ", "CE "},
["PHE"]={"CB ", "CG ", "CD1", "CD2", "CE1", "CE2", "CZ "},
["PRO"]={"CB ", "CG ", "CD "},
["SER"]={"CB ", "OG "},
["THR"]={"CB ", "OG1", "CG2"},
["TRP"]={"CB ", "CG ", "CD1", "CD2", "NE1", "CE2", "CE3", "CZ2", "CZ3", "CH2"},
["TYR"]={"CB ", "CG ", "CD1", "CD2", "CE1", "CE2", "CZ ", "OH "},
["VAL"]={"CB ", "CG1", "CG2"}
}


function mutate()

aa_chosen=document:getElementById("mutate").value
--print("New amino acid chosen:",aa_chosen)
--for k,l in ipairs(sidechain_atoms[aa_chosen]) do print(l) end
residue=protein[startpos]
resnam=string.sub(residue["CA "].res,1,3)
--print("Old amino acid:",resnam)
new_residue={}
for k,l in pairs(residue) do if k=="N  " or k=="CA " or k=="C  " or k=="O  " or (k=="CB " and aa_chosen~="GLY") then new_residue[k]=l end end  --copy N, CA, C, O and CB to new residue; exclude CB with mutations to GLY
--for m,n in pairs(new_residue) do print(m, "copied to new residue") end
aa_chosen=document:getElementById("mutate").value
--print("New amino acid chosen:",aa_chosen)
--for k,l in ipairs(sidechain_atoms[aa_chosen]) do print(l) end

car=residue["C  "]
cal=residue["CA "]
dx=cal.x-car.x
dy=cal.y-car.y
dz=cal.z-car.z

if resnam~="GLY" then cb=residue["CB "] else cb=residue["CA "] end

for k,l in ipairs(sidechain_atoms[aa_chosen]) do
    --print(l)
    if l~="CB " then
       newat={["x"]=cb.x+(k-1)*dx/1.5, ["y"]=cb.y+(k-1)*dy/1.5, ["z"]=cb.z+(k-1)*dz/1.5}
    new_residue[l]=newat
    elseif resnam=="GLY" then
    nit=residue["N  "]
    meanx=(nit.x+car.x)/2  meany=(nit.y+car.y)/2  meanz=(nit.z+car.z)/2  --adds CB to a GLY in plane with CA on opposite side of N and C (tetrig and checkChir should sort it).
    dmx=cal.x-meanx  dmy=cal.y-meany  dmz=cal.z-meanz
    newat={["x"]=cb.x+dmx*2, ["y"]=cb.y+dmy*2, ["z"]=cb.z+dmz*2}
    new_residue[l]=newat
    end
end

for j,k in pairs(new_residue) do
    --print("Atom:",j)
    k["res"]=aa_chosen
    k["num"]=cal.num
    k["chain"]=cal.chain
    k["nsymop"]={0}
    k["postrot_shift"]={0,0,0}
    --for l,m in pairs(k) do print(l,m) end
end

--for i,j in pairs(new_residue) do print(i) for k,l in pairs(j) do print(k,l) end end

residue=new_residue
resnam=string.sub(residue["CA "].res,1,3)
ncyc=0   --normally needed
sorted=false
repeat mmm() until sorted or ncyc>600
print("Residue:",resnam,startpos,"max delta:",maxdiff,"after", ncyc, "cycles.")
protein[startpos]=residue
document:getElementById("mutate").selectedIndex=0
startmove()
rerun() --needed normally

end



sidechbutton=document:getElementById("viewsidechain")
sidechbutton:addEventListener("click", function() if viewSideCh==false then sidechbutton.innerHTML="Main chain" viewSideCh=true run() else viewSideCh=false sidechbutton.innerHTML="Side chain" run() end end)

stepbutton=document:getElementById("step")
stepbutton:addEventListener("click", function() step(1) end)

stepbutton=document:getElementById("back")
stepbutton:addEventListener("click", function() step(-1) end)

stopbutton=document:getElementById("stop")
stopbutton:addEventListener("click", stop)

mutatebutton=document:getElementById("mutate")
mutatebutton:addEventListener("change", mutate)

movebutton=document:getElementById("move")
movebutton:addEventListener("click", startmove)

savebutton=document:getElementById("tidy")
savebutton:addEventListener("click", tidy)

undobutton=document:getElementById("undo")
undobutton:addEventListener("click", undo)

savebutton=document:getElementById("save")
savebutton:addEventListener("click", save)

pdbinput=document:getElementById("pdbInput")
pdbinput:addEventListener("change", function() readform() readfile() end)

mapinput=document:getElementById("mapInput")
mapinput:addEventListener("change", function() readform() readmapfile() end)

posneg=document:getElementById("densign")
posneg:addEventListener("click", swapdensign)

symtick=document:getElementById("symmetry")  --add callbacks to the text- and tick-boxes. 
symtick:addEventListener("change", run)

position=document:getElementById("position")
position:addEventListener("change", run)

mapcon=document:getElementById("contour")
mapcon:addEventListener("change", function() if minimap~=nil then if #(minimap)>0 then run() end end end)

canvas:addEventListener("mousedown", function(_,event) mouseDown(event) end)
canvas:addEventListener("mousemove", function(_,event) mouseMove(event) end)
canvas:addEventListener("mouseup",   function(_,event) mouseUp(event)   end)

canvas:addEventListener("touchstart", function(_,event) touchStart(event) end) -- add ",false" after "end" to stop browser reading swipe events?
canvas:addEventListener("touchmove", function(_,event) touchMove(event) end) -- ditto
canvas:addEventListener("touchend",  function(_,event) touchEnd(event)   end) -- ditto



-- Functions for drawing symmetry mates


function get_distance(x1,y1,z1,x2,y2,z2)
dx=x2-x1
dy=y2-y1
dz=z2-z1
d1=(a*dx)^2+(b*dy)^2+(c*dz)^2
d2=2*dx*dy*a*b*cos_gamma
d3=2*dy*dz*b*c*cos_alpha
d4=2*dz*dx*c*a*cos_beta
d=math.sqrt(d1+d2+d3+d4)
end


function get_centroid()
if readcell_ok and nCA~=0 then
    centx,centy,centz=xtot/nCA,ytot/nCA,ztot/nCA  --orthogonal
    cent_coord=window:Array(centx,centy,centz)
    fract_cent_coord=mathjs:multiply(Linv,cent_coord)  --fractional
    fcentx,fcenty,fcentz=fract_cent_coord:subset(el1),fract_cent_coord:subset(el2),fract_cent_coord:subset(el3)
    originx,originy,originz=math.floor(fcentx),math.floor(fcenty),math.floor(fcentz)
    oricoord=window:Array(originx,originy,originz)
    centroid=window:Array(fcentx,fcenty,fcentz)
    readlog="No. CA atoms read: "..nCA.."Orthogonal coordinates\n"
    readlog=readlog.."Min x,y,z:"..string.format("%12.7f %12.7f %12.7f\n",Xmin,Ymin,Zmin)
    readlog=readlog.."Max x,y,z:"..string.format("%12.7f %12.7f %12.7f\n",Xmax,Ymax,Zmax)
    readlog=readlog.."Centroid: "..string.format("%12.7f %12.7f %12.7f\n",centx,centy,centz)
    print(readlog)
    fill_symm_box()
elseif readcell_ok then print("No amino acids found.")
else print("No unit cell and symmetry found i.e. the CRYST1 line is missing.") end
end


function fill_symm_box()
if sg~="P 1       " then notP1=true else notP1=false end
rotbox,symm_box,tranlist={},{},{}
wideborder=20 border=boxrad  --border on orthogonal axes in Angstroms
nshift=0 chainout=0
steps=1 --tested fine with steps=1 for all 67 space groups below

for ashift=-steps,steps do  for bshift=-steps,steps do  for cshift=-steps,steps do table.insert(tranlist,{ashift,bshift,cshift}) end end end
total_steps=#(tranlist)
--print("Total number of steps:",total_steps)

if notP1 then 
	rotations()
else
	translations(protein)  --with P1 just do translations
end
end


function rotations()  --apply SG rotation and translation, get new centroid and new origin.
numrot=#(sg_rotations)
print("Number of rotations:",numrot)
window:setTimeout(function() document:getElementById("info").innerHTML="Please wait <progress value=\"0\" max=\"100\" id=\"progress\" style=\"width:100px;\"></progress>" end, 0)
rotations2(1)
end

function rotations2(k)  --rewritten as a recursion to get progress bar to work
rotop=sg_rotations[k]
rotrow1=window:Array(rotop[1][1],rotop[1][2],rotop[1][3])
rotrow2=window:Array(rotop[2][1],rotop[2][2],rotop[2][3])
rotrow3=window:Array(rotop[3][1],rotop[3][2],rotop[3][3])
rotrows=window:Array(rotrow1, rotrow2, rotrow3)
rotat=js.new(mathjs.matrix,rotrows)
trans=sg_translations[k]
trans_array=window:Array(trans[1],trans[2],trans[3])
orthrot=mathjs:chain(L):multiply(rotat):multiply(Linv):done()  --order is right
orthtrans=mathjs:multiply(L,trans_array)
newcent=mathjs:chain(rotat):multiply(centroid):add(trans_array):done()  --new centroid fractional coords
newxc,newyc,newzc=newcent:subset(el1),newcent:subset(el2),newcent:subset(el3)
newxo,newyo,newzo=math.floor(newxc),math.floor(newyc),math.floor(newzc)  --new origin
xshift=originx-newxo  yshift=originy-newyo  zshift=originz-newzo  --shift to move new molecule close to original
newnewxc=newxc+xshift  newnewyc=newyc+yshift  newnewzc=newzc+zshift  --new fractional centroid shifted close to original
dmin=math.huge

for oa=-2,2 do  --125 unit cell translations to get the "optimum" position (on a,b,c) of the symmetry mate (i.e. centroid is closest to the original)
	for ob=-2,2 do
		for oc=-2,2 do
			testxc=newnewxc+oa  testyc=newnewyc+ob  testzc=newnewzc+oc
			get_distance(testxc,testyc,testzc,fcentx,fcenty,fcentz)
			if d<dmin then dmin=d postrot_shift={xshift+oa,yshift+ob,zshift+oc} end  --postrot_shift is always fractional
		end
	end
end

fract_postrot_shift=window:Array(postrot_shift[1],postrot_shift[2],postrot_shift[3])  --need to orthogonalise the shifts
orth_post_rot_shift=mathjs:multiply(L,fract_postrot_shift)
postRotShiftx=orth_post_rot_shift:subset(el1)  postRotShifty=orth_post_rot_shift:subset(el2)  postRotShiftz=orth_post_rot_shift:subset(el3)
print("Orthogonal shifts needed after rotation:",k,"are:",postRotShiftx,postRotShifty,postRotShiftz)

for i,res in ipairs(protein) do  --&&
	residue={}
	for j,atom in pairs(res) do
		coord=window:Array(atom.x,atom.y,atom.z)
		res=atom.res num=atom.num chain=atom.chain  --copy out or protein table gets modified
		outcoord=mathjs:chain(orthrot):multiply(coord):add(orthtrans):add(orth_post_rot_shift):done()
		newx,newy,newz=outcoord:subset(el1),outcoord:subset(el2),outcoord:subset(el3)
		newatom={["res"]=res,["num"]=num,["chain"]=chain,["x"]=newx,["y"]=newy,["z"]=newz,["nsymop"]={k},["postrot_shift"]=postrot_shift}
		residue[j]=newatom
	end
	if next(residue)~=nil then table.insert(rotbox,residue) end
end  --&&

window:setTimeout(function()   --££
	prog=math.floor(k/numrot*100)
	progbar=document:getElementById("progress")
	progbar.value=tostring(prog)
	if k<numrot then  --##
		rotations2(k+1)  
	else
		--Note that no centroid juggling is done for the molecules generated by lattice centering (yet).
		if sg:sub(1,1)=="C" then centred_box={} censhift={f12,f12,0} centext={"+(&frac12;,&frac12;,0)"} addcentre(protein,1) addcentre(rotbox,1) for l,m in ipairs(centred_box) do table.insert(rotbox,m) end end
		if sg:sub(1,1)=="I" then centred_box={} censhift={f12,f12,f12} centext={"+(&frac12;,&frac12;,&frac12;)"} addcentre(protein,1) addcentre(rotbox,1) for l,m in ipairs(centred_box) do table.insert(rotbox,m) end end
		if sg:sub(1,1)=="F" then centred_box={} fcen={{f12,f12,0},{0,f12,f12},{f12,0,f12}} centext={"+(&frac12;,&frac12;,0)","+(0,&frac12;,&frac12;)","+(&frac12;,0,&frac12;)"} for n,o in ipairs(fcen) do censhift=o addcentre(protein,n) addcentre(rotbox,n) end for l,m in ipairs(centred_box) do table.insert(rotbox,m) end end
		translations(protein)
		translations(rotbox)
	end  --##
end, 0)  --££

end


function addcentre(thing,ncen)
centring_shift=window:Array(censhift[1],censhift[2],censhift[3])  --orthogonalise the centring shifts
orth_cen_shift=mathjs:multiply(L,centring_shift)
orth_cen_shift_x,orth_cen_shift_y,orth_cen_shift_z=orth_cen_shift:subset(el1),orth_cen_shift:subset(el2),orth_cen_shift:subset(el3)
for i,res in ipairs(thing) do
    residue={}
    for j,atom in pairs(res) do
        newx,newy,newz=atom.x+orth_cen_shift_x,atom.y+orth_cen_shift_y,atom.z+orth_cen_shift_z
        res=atom.res num=atom.num chain=atom.chain nsymop=atom.nsymop[1] postrot_shift=atom.postrot_shift --print(nsymop) --copy out or protein table gets modified
        newatom={["res"]=res,["num"]=num,["chain"]=chain,["x"]=newx,["y"]=newy,["z"]=newz,["nsymop"]={nsymop,ncen},["postrot_shift"]=postrot_shift}
        residue[j]=newatom
    end
    if next(residue)~=nil then table.insert(centred_box,residue) end
end
end


function translations(thing)
window:setTimeout(function() document:getElementById("info").innerHTML="Please wait <progress value=\"0\" max=\"100\" id=\"progress\" style=\"width:100px;\"></progress>" end, 0)
translations2(thing,1)
end

function translations2(thing,num_tran)  --rewritten as a recursion to get out of progress bar hell
ashift,bshift,cshift=tranlist[num_tran][1],tranlist[num_tran][2],tranlist[num_tran][3]
--print(ashift,bshift,cshift)
abc_shift=window:Array(ashift,bshift,cshift)
orth_abc_shift=mathjs:multiply(L,abc_shift)
orth_abc_shift_x,orth_abc_shift_y,orth_abc_shift_z=orth_abc_shift:subset(el1),orth_abc_shift:subset(el2),orth_abc_shift:subset(el3)
nshift=nshift+1
--print("nshift",nshift)
if thing==protein and ashift==0 and bshift==0 and cshift==0 then do end else      --#
	for i,res in ipairs(thing) do   --$
		residue={}
		_,fastatom=next(res)  --first (i.e. random) atom in residue must be close to the border (within wideborder) before doing atom-by-atom scan
		newx=fastatom.x+orth_abc_shift_x  newy=fastatom.y+orth_abc_shift_y newz=fastatom.z+orth_abc_shift_z
		if Xmin-wideborder<=newx and newx<=Xmax+wideborder and Ymin-wideborder<=newy and newy<=Ymax+wideborder and Zmin-wideborder<=newz and newz<=Zmax+wideborder then   --!"
			for j,atom in pairs(res) do   --£
				newx=atom.x+orth_abc_shift_x  newy=atom.y+orth_abc_shift_y  newz=atom.z+orth_abc_shift_z
				if Xmin-border<=newx and newx<=Xmax+border and Ymin-border<=newy and newy<=Ymax+border and Zmin-border<=newz and newz<=Zmax+border then
					res=atom.res num=atom.num chain=atom.chain nsymop=atom.nsymop postrot_shift=atom.postrot_shift  --copy out or protein table gets modified
					new_postrot_shift={atom.postrot_shift[1]+ashift, atom.postrot_shift[2]+bshift, atom.postrot_shift[3]+cshift}
					--print("Inside box:",atom.res,atom.num,ashift,bshift,cshift)
					--if oldnshift~=nshift then oldnshift=nshift chainout=chainout+1 print("chainout",chainout) end
					newatom={["res"]=res,["num"]=num,["chain"]=chain,["x"]=newx,["y"]=newy,["z"]=newz,["nsymop"]=nsymop,["postrot_shift"]=new_postrot_shift}
					residue[j]=newatom
				end
			end  --£
			if next(residue)~=nil then table.insert(symm_box,residue) end
		end --!"
	end  --$
end --#

prog=math.floor(num_tran/total_steps*100)

if num_tran<total_steps then
	window:setTimeout(function()
	progbar=document:getElementById("progress")
	progbar.value=tostring(prog)
	translations2(thing,num_tran+1) 
	end, 0)
elseif thing==rotbox or notP1==false then
	symm_box_filled=true
	run()  --bit ugly to run all of run() but seems necessary when things running asynchronously - more progress bar hell 
end

end


f16=1/6 f14=1/4 f13=1/3 f12=1/2 f23=2/3 f34=3/4 f56=5/6
frac_lookup={[f16]="&frac16;",[f14]="&frac14;",[f13]="&frac13;",[f12]="&frac12;",[f23]="&frac23;",[f34]="&frac34;",[f56]="&frac56;"}

symops={
["P 1 2 1   "]={rotas={{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,0,0}}},
["P 1 21 1  "]={rotas={{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,f12,0}}},
["C 1 2 1   "]={rotas={{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,0,0}}},
["I 1 2 1   "]={rotas={{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,0,0}}},

["P 21 21 21"]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{f12,f12,0},{0,f12,f12},{f12,0,f12}}},
["P 21 21 2 "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{f12,f12,0},{f12,f12,0},{0,0,0}}},
["P 2 21 21 "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{0,0,0},{0,f12,f12},{0,f12,f12}}},
["P 2 2 2   "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{0,0,0},{0,0,0},{0,0,0}}},
["P 2 2 21  "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{0,0,0},{0,0,f12},{0,0,f12}}},
["P 21 2 21 "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{f12,0,f12},{0,0,0},{f12,0,f12}}},

["C 2 2 21  "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{0,0,0},{0,0,f12},{0,0,f12}}},
["C 2 2 2   "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{0,0,0},{0,0,0},{0,0,0}}},
["I 2 2 2   "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{0,0,0},{0,0,0},{0,0,0}}},
["I 21 21 21"]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{f12,f12,0},{0,f12,f12},{f12,0,f12}}},
["F 2 2 2   "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}}},trans={{0,0,0},{0,0,0},{0,0,0}}},

["P 4       "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,0},{0,0,0},{0,0,0}}},
["P 41      "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,f14},{0,0,f12},{0,0,f34}}},
["P 42      "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,f12},{0,0,0},{0,0,f12}}},
["P 43      "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,f34},{0,0,f12},{0,0,f14}}},
["P 4 2 2   "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},
["P 41 2 2  "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,0,f14},{0,0,f12},{0,0,f34},{0,0,f34},{0,0,f12},{0,0,f14},{0,0,0}}},
["P 42 2 2  "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,0,f12},{0,0,0},{0,0,f12},{0,0,f12},{0,0,0},{0,0,f12},{0,0,0}}},
["P 43 2 2  "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,0,f34},{0,0,f12},{0,0,f14},{0,0,f14},{0,0,f12},{0,0,f34},{0,0,0}}},
["P 4 21 2  "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{f12,f12,0},{0,0,0},{f12,f12,0},{0,0,0},{f12,f12,0},{0,0,0},{f12,f12,0}}},
["P 41 21 2 "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{f12,f12,f14},{0,0,f12},{f12,f12,f34},{0,0,0},{f12,f12,f34},{0,0,f12},{f12,f12,f14}}},
["P 42 21 2 "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{f12,f12,f12},{0,0,0},{f12,f12,f12},{0,0,0},{f12,f12,f12},{0,0,0},{f12,f12,f12}}},
["P 43 21 2 "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{f12,f12,f34},{0,0,f12},{f12,f12,f14},{0,0,0},{f12,f12,f14},{0,0,f12},{f12,f12,f34}}},
["I 4       "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,0},{0,0,0},{0,0,0}}},
["I 41      "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}}},trans={{0,f12,f14},{f12,f12,f12},{f12,0,f34}}},
["I 4 2 2   "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},
["I 41 2 2  "]={rotas={{{0,-1,0},{1,0,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}}},trans={{0,f12,f14},{f12,f12,f12},{f12,0,f34},{0,0,0},{f12,0,f34},{f12,f12,f12},{0,f12,f14}}},

["P 6       "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},
["P 61      "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}}},trans={{0,0,f16},{0,0,f13},{0,0,f12},{0,0,f23},{0,0,f56}}},
["P 62      "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}}},trans={{0,0,f13},{0,0,f23},{0,0,0},{0,0,f13},{0,0,f23}}},
["P 63      "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}}},trans={{0,0,f12},{0,0,0},{0,0,f12},{0,0,0},{0,0,f12}}},
["P 64      "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}}},trans={{0,0,f23},{0,0,f13},{0,0,0},{0,0,f23},{0,0,f13}}},
["P 65      "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}}},trans={{0,0,f56},{0,0,f23},{0,0,f12},{0,0,f13},{0,0,f16}}},

["P 6 2 2   "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},
["P 61 2 2  "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}}},trans={{0,0,f16},{0,0,f13},{0,0,f12},{0,0,f23},{0,0,f56},{0,0,f13},{0,0,f16},{0,0,0},{0,0,f56},{0,0,f23},{0,0,f12}}},
["P 62 2 2  "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}}},trans={{0,0,f13},{0,0,f23},{0,0,0},{0,0,f13},{0,0,f23},{0,0,f23},{0,0,f13},{0,0,0},{0,0,f23},{0,0,f13},{0,0,0}}},
["P 63 2 2  "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}}},trans={{0,0,f12},{0,0,0},{0,0,f12},{0,0,0},{0,0,f12},{0,0,0},{0,0,f12},{0,0,0},{0,0,f12},{0,0,0},{0,0,f12}}},
["P 64 2 2  "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}}},trans={{0,0,f23},{0,0,f13},{0,0,0},{0,0,f23},{0,0,f13},{0,0,f13},{0,0,f23},{0,0,0},{0,0,f13},{0,0,f23},{0,0,0}}},
["P 65 2 2  "]={rotas={{{1,-1,0},{1,0,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{-1,1,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}}},trans={{0,0,f56},{0,0,f23},{0,0,f12},{0,0,f13},{0,0,f16},{0,0,f23},{0,0,f56},{0,0,0},{0,0,f16},{0,0,f13},{0,0,f12}}},

["P 3       "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,0},{0,0,0}}},
["P 31      "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,f13},{0,0,f23}}},
["P 32      "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,f23},{0,0,f13}}},

["P 3 2 1   "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},
["P 31 2 1  "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}}},trans={{0,0,f13},{0,0,f23},{0,0,0},{0,0,f23},{0,0,f13}}},
["P 32 2 1  "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}}},trans={{0,0,f23},{0,0,f13},{0,0,0},{0,0,f13},{0,0,f23}}},

["P 3 1 2   "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},
["P 31 1 2  "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}}},trans={{0,0,f13},{0,0,f23},{0,0,f23},{0,0,f13},{0,0,0}}},
["P 32 1 2  "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{-1,1,0},{0,1,0},{0,0,-1}},{{1,0,0},{1,-1,0},{0,0,-1}}},trans={{0,0,f23},{0,0,f13},{0,0,f13},{0,0,f23},{0,0,0}}},

["H 3       "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}}, {{1,0,0},{0,1,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},  {{1,0,0},{0,1,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}}},trans={{0,0,0},{0,0,0},{f23,f13,f13},{f23,f13,f13},{f23,f13,f13}, {f13,f23,f23},{f13,f23,f23},{f13,f23,f23} }},
["H 3 2     "]={rotas={{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}},{{1,0,0},{0,1,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}},{{1,0,0},{0,1,0},{0,0,1}},{{0,-1,0},{1,-1,0},{0,0,1}},{{-1,1,0},{-1,0,0},{0,0,1}},{{0,1,0},{1,0,0},{0,0,-1}},{{1,-1,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{-1,1,0},{0,0,-1}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{f23,f13,f13},{f23,f13,f13},{f23,f13,f13},{f23,f13,f13},{f23,f13,f13},{f23,f13,f13},{f13,f23,f23},{f13,f23,f23},{f13,f23,f23},{f13,f23,f23},{f13,f23,f23},{f13,f23,f23}}},

["P 2 3     "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},
["P 21 3    "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}}},trans={{f12,f12,0},{0,f12,f12},{f12,0,f12},{0,0,0},{f12,0,f12},{f12,f12,0},{0,f12,f12},{0,0,0},{0,f12,f12},{f12,0,f12},{f12,f12,0}}},
["I 2 3     "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},
["I 21 3    "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}}},trans={{f12,f12,0},{0,f12,f12},{f12,0,f12},{0,0,0},{f12,0,f12},{f12,f12,0},{0,f12,f12},{0,0,0},{0,f12,f12},{f12,0,f12},{f12,f12,0}}},
["F 2 3     "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}}},trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},

["P 4 3 2   "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}},
      {{0,1,0},{1,0,0},{0,0,-1}},{{0,-1,0},{1,0,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{0,0,1},{0,1,0},{-1,0,0}},{{0,0,-1},{0,-1,0},{-1,0,0}},{{0,0,-1},{0,1,0},{1,0,0}},{{0,0,1},{0,-1,0},{1,0,0}},{{1,0,0},{0,0,1},{0,-1,0}},{{1,0,0},{0,0,-1},{0,1,0}},{{-1,0,0},{0,0,-1},{0,-1,0}},{{-1,0,0},{0,0,1},{0,1,0}}},
      trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},

["P 41 3 2  "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}},
      {{0,1,0},{1,0,0},{0,0,-1}},{{0,-1,0},{1,0,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{0,0,1},{0,1,0},{-1,0,0}},{{0,0,-1},{0,-1,0},{-1,0,0}},{{0,0,-1},{0,1,0},{1,0,0}},{{0,0,1},{0,-1,0},{1,0,0}},{{1,0,0},{0,0,1},{0,-1,0}},{{1,0,0},{0,0,-1},{0,1,0}},{{-1,0,0},{0,0,-1},{0,-1,0}},{{-1,0,0},{0,0,1},{0,1,0}}},
      trans={{f12,f12,0},{0,f12,f12},{f12,0,f12},{0,0,0},{f12,0,f12},{f12,f12,0},{0,f12,f12},{0,0,0},{0,f12,f12},{f12,0,f12},{f12,f12,0},{f34,f14,f14},{f14,f34,f14},{f14,f14,f34},{f34,f34,f34},{f34,f14,f14},{f34,f34,f34},{f14,f34,f14},{f14,f14,f34},{f34,f14,f14},{f14,f14,f34},{f34,f34,f34},{f14,f34,f14}}},

["P 42 3 2  "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}},
      {{0,1,0},{1,0,0},{0,0,-1}},{{0,-1,0},{1,0,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{0,0,1},{0,1,0},{-1,0,0}},{{0,0,-1},{0,-1,0},{-1,0,0}},{{0,0,-1},{0,1,0},{1,0,0}},{{0,0,1},{0,-1,0},{1,0,0}},{{1,0,0},{0,0,1},{0,-1,0}},{{1,0,0},{0,0,-1},{0,1,0}},{{-1,0,0},{0,0,-1},{0,-1,0}},{{-1,0,0},{0,0,1},{0,1,0}}},
      trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12},{f12,f12,f12}}},

["P 43 3 2  "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}},
      {{0,1,0},{1,0,0},{0,0,-1}},{{0,-1,0},{1,0,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{0,0,1},{0,1,0},{-1,0,0}},{{0,0,-1},{0,-1,0},{-1,0,0}},{{0,0,-1},{0,1,0},{1,0,0}},{{0,0,1},{0,-1,0},{1,0,0}},{{1,0,0},{0,0,1},{0,-1,0}},{{1,0,0},{0,0,-1},{0,1,0}},{{-1,0,0},{0,0,-1},{0,-1,0}},{{-1,0,0},{0,0,1},{0,1,0}}},
      trans={{f12,f12,0},{0,f12,f12},{f12,0,f12},{0,0,0},{f12,0,f12},{f12,f12,0},{0,f12,f12},{0,0,0},{0,f12,f12},{f12,0,f12},{f12,f12,0},{f14,f34,f34},{f34,f14,f34},{f34,f34,f14},{f14,f14,f14},{f14,f34,f34},{f14,f14,f14},{f34,f14,f34},{f34,f34,f14},{f14,f34,f34},{f34,f34,f14},{f14,f14,f14},{f34,f14,f34}}},

["I 4 3 2   "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}},
      {{0,1,0},{1,0,0},{0,0,-1}},{{0,-1,0},{1,0,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{0,0,1},{0,1,0},{-1,0,0}},{{0,0,-1},{0,-1,0},{-1,0,0}},{{0,0,-1},{0,1,0},{1,0,0}},{{0,0,1},{0,-1,0},{1,0,0}},{{1,0,0},{0,0,1},{0,-1,0}},{{1,0,0},{0,0,-1},{0,1,0}},{{-1,0,0},{0,0,-1},{0,-1,0}},{{-1,0,0},{0,0,1},{0,1,0}}},
      trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},

["I 41 3 2  "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}},
      {{0,1,0},{1,0,0},{0,0,-1}},{{0,-1,0},{1,0,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{0,0,1},{0,1,0},{-1,0,0}},{{0,0,-1},{0,-1,0},{-1,0,0}},{{0,0,-1},{0,1,0},{1,0,0}},{{0,0,1},{0,-1,0},{1,0,0}},{{1,0,0},{0,0,1},{0,-1,0}},{{1,0,0},{0,0,-1},{0,1,0}},{{-1,0,0},{0,0,-1},{0,-1,0}},{{-1,0,0},{0,0,1},{0,1,0}}},
      trans={{f12,f12,0},{0,f12,f12},{f12,0,f12},{0,0,0},{f12,0,f12},{f12,f12,0},{0,f12,f12},{0,0,0},{0,f12,f12},{f12,0,f12},{f12,f12,0},{f14,f34,f34},{f34,f14,f34},{f34,f34,f14},{f14,f14,f14},{f14,f34,f34},{f14,f14,f14},{f34,f14,f34},{f34,f34,f14},{f14,f34,f34},{f34,f34,f14},{f14,f14,f14},{f34,f14,f34}}},

["F 4 3 2   "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}},
      {{0,1,0},{1,0,0},{0,0,-1}},{{0,-1,0},{1,0,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{0,0,1},{0,1,0},{-1,0,0}},{{0,0,-1},{0,-1,0},{-1,0,0}},{{0,0,-1},{0,1,0},{1,0,0}},{{0,0,1},{0,-1,0},{1,0,0}},{{1,0,0},{0,0,1},{0,-1,0}},{{1,0,0},{0,0,-1},{0,1,0}},{{-1,0,0},{0,0,-1},{0,-1,0}},{{-1,0,0},{0,0,1},{0,1,0}}},
      trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0}}},

["F 41 3 2  "]={rotas={{{1,0,0},{0,-1,0},{0,0,-1}},{{-1,0,0},{0,1,0},{0,0,-1}},{{-1,0,0},{0,-1,0},{0,0,1}},{{0,1,0},{0,0,1},{1,0,0}},{{0,-1,0},{0,0,-1},{1,0,0}},{{0,1,0},{0,0,-1},{-1,0,0}},{{0,-1,0},{0,0,1},{-1,0,0}},{{0,0,1},{1,0,0},{0,1,0}},{{0,0,-1},{1,0,0},{0,-1,0}},{{0,0,-1},{-1,0,0},{0,1,0}},{{0,0,1},{-1,0,0},{0,-1,0}},
      {{0,1,0},{1,0,0},{0,0,-1}},{{0,-1,0},{1,0,0},{0,0,1}},{{0,1,0},{-1,0,0},{0,0,1}},{{0,-1,0},{-1,0,0},{0,0,-1}},{{0,0,1},{0,1,0},{-1,0,0}},{{0,0,-1},{0,-1,0},{-1,0,0}},{{0,0,-1},{0,1,0},{1,0,0}},{{0,0,1},{0,-1,0},{1,0,0}},{{1,0,0},{0,0,1},{0,-1,0}},{{1,0,0},{0,0,-1},{0,1,0}},{{-1,0,0},{0,0,-1},{0,-1,0}},{{-1,0,0},{0,0,1},{0,1,0}}},
      trans={{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{0,0,0},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14},{f14,f14,f14}}},
}

</script>
</body>
</html>

minimaphelp.html

htmla year ago

Help file for minimapai.de

<html>
<head><meta name="viewport" content="width=device-width, initial-scale=1.0">
<title>mini map aide help</title>
<style>
body {margin:auto; font-size:16px; font-family:arial; color:dodgerblue; background-color:#0C0A0A;}
p {margin:10px;}
</style>
</head>
<body>
<p>
<a href="index.html"><b style="color:#FFD633;">Mini Map Aide</b></a> allows the fit of a protein structure, in PDB format, to its electron density map, in CCP4 format, to be assessed using a web browser on mobile devices. It seems to work with nucleotides and small ligands, too.
<p>
It is meant to be a completely minimal map-checking aid that works on a residue-by-residue basis and could be used for checking a structure in the final refinement rounds whenever the user has no computer access. The bad regions could be recorded, e.g. by taking screen-shots, for later attention with <a href="https://www2.mrc-lmb.cam.ac.uk/personal/pemsley/coot/">Coot</a>. Despite its limitations, the aim is for it to be useful in research or as a teaching aid.
<p>
If you have a PDB file and the corresponding map, the display calculations are all done by a client-side script in the browser and neither file is uploaded to the server.
<p>
The map must be in CCP4 format and <u>must</u> either cover the protein molecule or fill an entire unit cell. On a computer with CCP4 installed, this can be achieved most easily from refinement MTZ files with the following command: <p style="color:#FFD633;">gemmi sf2map -s3 --mapmask=my.pdb my.mtz my.map</p><p> or requires running MAPMASK in the old CCP4i GUI or using the Phenix Maps tab. With gemmi the -s number can be increased for finer grid sampling. Use the <span style="color:#FFD633;">-d</span> option to make a difference map from a refmac MTZ file. Note, there are some <a href="https://gemmi.readthedocs.io/en/latest/utils.html#sf2map">constraints</a> on the expected MTZ column labels.
<p>
CCP4 maps for most proteins can be obtained from <a href="https://www.ebi.ac.uk/pdbe/">EBI PDBe</a> in "Downloads" as "EDS map" and "EDS difference map."
<p>
<span style="color:#FFD633;">Mini Map Aide can also calculate the map from an MTZ file.</span> For this just use the map file chooser to pick the MTZ file instead. The PDB and reflection files will then be uploaded to the server to make the map file which you can then download and save. The map must then be read back into the program in order for it to be displayed. The files on the server will be deleted a few minutes after they are created.
<p>
The contour level can be adjusted and is in units of the map sigma or RMSD. If a negative value is used, the contour lines become red.
<p>
More help is available in the tooltips which appear on a PC browser when you hover the mouse over any of the text.
<p>
If it is too slow, try increasing the contour level to reduce the size of the mesh.
<p>
Things it does do (I hope):
<ul>
		<li>Makes maps from MTZ files using Gemmi by Marcin Wojdyr (CCP4/Global Phasing).</li>
		<li>Handles mono- and triclinic maps, where <b>b</b> is along Z (ncode 3).</li>
		<li>Shows negative contours in red, so one can look at difference maps. </li>
		<li>Draws multiple conformations and disulphides (the latter in pale green).</li>
		<li>Displays <i>small</i> EM maps from <a href="https://www.ebi.ac.uk/emdb">EMDB</a>.</li>
		<li>Allows small changes to the amino acids to be made.</li>
		<li>Basic regularisation of amino acid geometry, including C&alpha; and C&beta; chirality.</li>
                <li>Allows amino acid mutations to be made.</li>
                <li>Saving a new PDB file.</li>
		<li>Draws symmetry mates.</li>
</ul>
<p>
Things it <b>cannot</b> do:
<ul>
       <li>Display more than one map at a time.</li>
</ul>
<p>
<b>Input fields</b>
<p>
<b>Go to:</b> The current position in the sequence. As soon as the PDB file is read, this will allow the user to choose which chain and residue number to centre on. Click redraw to allow any changes to take effect. If a large positive number is entered for a given chain, e.g. A 99999, the last residue in that chain will be displayed, which might be a a water molecule. Likewise, a negative number e.g. A -9999, would cause the program to jump back to the first residue in that chain, e.g. A -7. If the residue chosen is actually missing, e.g. in a disordered loop, the program will centre on the next available amino acid. Each residue will be displayed with its CA atom at the centre of rotation and its label at the top of the graphics window. If the residue possesses a CB atom, the electron density map will be centred on this, instead of the CA, to aid viewing the side chain fit. With nucleotides, the default centering is on the P atom.
<p>
<b>Contour level:</b> The contour level in sigma or rms units. Change sign with +/- button.
<p>
<b>PDB:</b> and <b>Map:</b> are file choosers for the coordinates and map or MTZ files.
<p>
<b>Side chain:</b> Click this button to centre on the outer parts of amino acid side chains or on nucleotide bases. Clicking toggles between centering on the side chain and the main chain C&alpha;.
<p>
<b>Roll / Scroll:</b> This switches between animation and scrolling of the molecule.
<p>
<b>Last residue:</b> Move to preceding residue in the sequence.
<p>
<b>Next residue:</b> Move to next residue in the sequence.
<p>
<b>Symmetry:</b> Ticking this checkbox causes the nearby parts of symmetry-related molecules to be shown as green sticks. It can take a while to calculate all the symmetry mates the first time this option is selected, so a progress bar or two might appear and patience is needed, but it will be quicker thereafter. It works for nearly all of the 65 enantiomorphic space groups and a few other non-standard ones that occur in the PDB.
<p>
<b>Mutate:</b> This drop-down menu allows the amino acid at the centre of the display to be mutated. The side chain of the new amino acid will probably require rebuilding by hand using the Move and Tidy options. If the outcome is satisfactory, the Save option will be needed to write the coordinates to a new local PDB file.
<p>
<b>Move:</b> This allows atoms in the residue at the centre of the display to be moved. The moveable atoms will appear green and can be clicked and dragged. The atom selected for movement is enlarged momentarily to indicate which one is actually moving. The user might find that just moving some of the side chain atoms in the right direction followed by use of the Tidy function may give an improved fit to the map. Once the Move button is clicked, the following 3 buttons become active.
<p>
<b>Tidy:</b> This should improve the geometry of an amino acid after any of the atoms have been moved and redraw the bonds. It is rather approximate although the final bond lengths and interatomic distances defining the bond angles should be correct to within about 0.1 <span>&#8491;</span>. It does not use the electron density map at all but it may be repeated as many times as required since each run can give a slightly different, but not necessarily better, fit. Note that in order to retain reasonable main chain geometry it does not move the backbone N and C atoms, so these are probably best not moved too far by the user.
<p>
<b>Undo:</b> This should completely undo the move and tidy operations. Note that unlike, Word, Coot, <i>etc</i>, the undo operations are not chained <i>i.e.</i> only the current round of changes made to the residue at the centre of the screen can be undone.
<p>
<b>Save:</b> This allows the changes to be saved in a new PDB file. It is recommended to do this at least once before closing the browser tab if the results are wanted for future use! Note that if no changes have been made to the residue currently at the screen centre, making the Save button active would require clicking the Move button first, although it would not be necessary to move any of the atoms.
<p>
The script (which can be downloaded within the html of the main page) is written in <a href="https://www.lua.org/">Lua</a> and uses <a href="http://fengari.io/">Fengari</a> to communicate with <a href="https://mathjs.org/">Math.js</a> and <a href="https://threejs.org/">Three.js</a>. The contouring algorithm (surfacenets.js), by Mikola Lysenko, is described <a href="https://www.merl.com/publications/docs/TR99-24.pdf">here</a> and was obtained from <a href="https://github.com/mikolalysenko/isosurface">here</a>. Help from the <a href="https://www.reddit.com/r/lua/">Lua reddit</a> group was invaluable.
<p>
I am very grateful to Marcin Wodjyr (CCP4/Global Phasing) for allowing me to use <a href="https://gemmi.readthedocs.io/en/latest/utils.html#sf2map">gemmi sf2map</a>.
<p>
Please report any problems or corrections to jbcooper_at_fastmail_dot_net who will be glad to help.

crontab file for minimapai.de

plain_text2 years ago

*/5  * * * * find /var/www/html/temp -type f -mmin +5 -exec rm "{}" \;
0 0 * * * certbot --webroot renew

Server script (run_gemmi.py) for minimapai.de

python4 years ago

#!/usr/bin/python3
import cgi
import os
import cgitb
import subprocess
print("Content-type:text/html\r\n\r\n") # needed
cgitb.enable()
form = cgi.FieldStorage()
pdbfile = form["pdbfile"]
pdbname = os.path.basename(pdbfile.filename)
open("/var/www/html/temp/" + pdbname, "wb").write(pdbfile.file.read())
#message = "<p>The file " + pdbname + " was uploaded successfully</p>"
#print(message)
mtzfile = form["mtzfile"]
mtzname = os.path.basename(mtzfile.filename)
open("/var/www/html/temp/" + mtzname, "wb").write(mtzfile.file.read())
#message = "<p>The file " + mtzname + " was uploaded successfully</p>"
#print(message)
diffneeded = form.getvalue("diffneeded")
#print(diffneeded)
mapname=mtzname.split(".",1)[0]
if diffneeded=="true": gemmiproc=["gemmi","sf2map","-d","-s3","--mapmask=../temp/"+pdbname,"../temp/"+mtzname,"../temp/"+mapname+"_diff.map"]
else: gemmiproc=["gemmi","sf2map","-s3","--mapmask=../temp/"+pdbname,"../temp/"+mtzname,"../temp/"+mapname+".map"]
#print("<p>",gemmiproc,"</p>")
gemmirun=subprocess.run(gemmiproc)

Total 4 snippets
1
5 / page