dm-script - 如何移动 SI 图像中的列

Question

为了在SI图像中进行漂移校正，如下图所示:

我写代码:

number max_shift=5
image src := GetFrontImage()
number sx, sy, sz
src.Get3DSize(sx, sy, sz)
result("sx: "+sx+"\n")
result("sy: "+sy+"\n")
result("sz: "+sz+"\n")

// assume a random shift in x
image shift := IntegerImage("xcorrd",4,0, 1, sy)
shift = max_shift*Random() 

// make a coordinate table
image col := IntegerImage("col",4,0, sx, sy)
image row := IntegerImage("row",4,0, sx, sy)
image plane := IntegerImage("plane",4,0, sx, sy)
col = icol
row = irow
plane = iplane

// to expand the shift as the same size with source image
image ones := IntegerImage("ones",4,0, sx, sy)
ones = 1

// create a random column shift of the source SI image
for(number i=0; i<sy; i++) {
    col[i,0,i+1,sx] = col[i,0,i+1,sx]+shift.GetPixel(0,i)*ones[i,0,i+1,sx]
};

// drift corrected
image im := RealImage("test si", 4, sx+max_shift, sy, sz)
im=0
im[col, row, plane] = src[icol,irow,iplane]

im.ImageGetTagGroup().TagGroupCopyTagsFrom(src.ImageGetTagGroup())
im.ImageCopyCalibrationFrom(src)
im.SetName(src.GetName()+"-drift corrected")
im.showimage()

可以校正图像，但是无法将光谱传输到校正后的 SI，如图所示:

我只是想知道我的脚本有什么问题。

提前谢谢你。

Answer 1

im[col, row, plane] = src[icol,irow,iplane]

内部变量 icol、irow、iplane 将由行中唯一的固定大小图像表达式求值。在你的情况下，行和平面(所有大小相同)

但是，它们都是二维的，所以内部发生的是您迭代 X 和 Y，然后写入值:

im[ col(x,y), row(x,y), plane(x,y) ] = src[x,y,0]//遍历所有 x/y

正如 Don I 在评论中提到的那样，您可能希望迭代 z 维度。

或者，您可以在脚本中制作所有大小为 (sx,sy,sz) 的图像。这适用于表达式，但效率极低。

一般来说，这里最好的解决方案是根本不使用 icol、irow、iplane，而是使用 Slice 命令。见this answer :

---

我可能会像下面这样为 SI 编写一个逐行的 x 位移:该脚本利用了这样一个事实，即可以在 x 方向上移动整个“ block ”(X x 1 x Z)，在 y 上迭代。

number sx = 256
number sy = 256
number sz = 100
image testSI := realImage("SI",4,sx,sy,sz)
testSI = sin(itheta/(idepth-iplane)*idepth) + (iplane % (icol+1))/idepth
testSI.ShowImage()

image xdrift := RealImage("XDrift per line",4,sy)
xdrift = trunc(random()*5 + 20*sin(icol/iwidth*3*PI()))
xdrift.ShowImage()

// Apply linewise Drift to SI, assuming xDrift holds this data
xDrift -= min(xDrift)   // ensure only positive shifts
image outSI := realImage("SI shifted",4,sx+max(xDrift),sy,sz)
outSI.ShowImage()

for( number y=0; y<sy; y++ ){
    number yShift = sum(xDrift[y,0])
    outSI.slice2( yShift,y,0, 0,sx,1, 2,sz,1 ) = testSI.slice2(0,y,0,0,sx,1,2,sz,1)
}

---

下面的脚本执行“逐个平面”的迭代，但对平面移位没有限制。事实上，这里每个像素都有一个指定的 XY 位移。

请注意，如果您想使用值的双线性插值，则可以使用 warp(source, xexpr, yexpr ) 而不是二维寻址 source[ xexpr, yexpr ] (和有效范围外的 0 截断)。

number sx = 256
number sy = 256
number sz = 100
image testSI := realImage("SI",4,sx,sy,sz)
testSI = sin(itheta/(idepth-iplane)*idepth) + (iplane % (icol+1))/idepth
testSI.ShowImage()

image xdrift := RealImage("pixelwise XDrift",4,sx,sy)
xdrift = irow%10*random() + 20*cos(irow/iheight*5*PI())
xdrift.ShowImage()

image ydrift := RealImage("pixelwise yDrift",4,sx,sy)
ydrift = 10*abs(cos(icol/iwidth* (irow/iheight) * 10 * PI())) + irow/iheight * 10 
ydrift.ShowImage()

// Apply pixelwise Drift to SI
xDrift -= min(xDrift)   // ensure only positive shifts
yDrift -= min(yDrift)   // ensure only positive shifts
number mx = max(xDrift)
number my = max(yDrift)
image outSI := realImage("SI shifted",4,sx+mx,sy+my,sz)
outSI.ShowImage()
for( number z=0;z<sz;z++){
    image outPlane := outSI.Slice2( 0,0,z, 0,sx+mx,1,1,sy+my,1) 
    image srcPlane := testSI.Slice2( 0,0,z, 0,sx,1,1,sy,1)
    outPlane = srcPlane[ icol + xdrift[icol,irow] - mx, irow + ydrift[icol,irow] - my ]
    // outPlane = srcPlane.warp( icol + xdrift[icol,irow] - mx, irow + ydrift[icol,irow] - my)
}