how to interpolate missing data within a 3D cube of data, or for a 2D image with missing slices
nteutsch
Wed, 10/25/2017 - 03:22 am
I am working with a set number 2D images each taken at a different z values (-7 up to +8, steps of 1). I have been able to combine these images into a 3D cube with the z axis being evenly spaced and then used gizmo to view the "surfaces" of the box, no problem.
My problem comes with my second data set, which has missing voltages (was only taken at -7,-5,-3,-1,0,2,4,6,8). I have made this into a 3D cube again, with even spacing again, and have essentially filled in the missing steps with 2D images of NaN values.
My question is: how can I interpolate to fill in this missing data slices?
Currently I am taking a slice in the y,z plane for each given x value, to create a 2D image and am then trying to use Interp2D to fill in the NaN slices but this isn't working (in the second half of the code provided below). This is the image I have attached. How can I interpolate this?
I would then re-populate the 3D cube with the interpolated data.
Thanks in advance for any help!
**********************
Function make_3D_cube()
//make waves to hold the information for the slices and correponding volatges
string slices="SMP1212_028_kx;temp_neg6;SMP1212_042_kx;temp_neg4;SMP1212_041_kx;temp_neg2;SMP1212_039_kx;SMP1212_037_kx;temp_1;SMP1212_036_kx;temp_3;SMP1212_035_kx;temp_5;SMP1212_034_kx;temp_7;SMP1212_033_kx"
Make/O voltages={-7, -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8}
String k_scale_str=stringfromlist(0,slices,";")
wave k_scale = $k_scale_str
//make temp waves to fill in missing voltages with NaN waves
duplicate/O k_scale, temp_7
temp_7=nan
duplicate/O temp_7, temp_neg6, temp_neg4, temp_neg2, temp_1, temp_3, temp_5
//make 3D cube to poulate
Make/O/N=(DimSize(k_scale,1), DimSize(k_scale,0), DimSize(voltages,0)) cube_3D
SetScale/P x, DimOffset(k_scale,1), DimDelta(k_scale,1), "Å^-1", cube_3D
SetScale/P y, DimOffset(k_scale,0), DimDelta(k_scale,0), "eV", cube_3D
SetScale/I z, voltages[0], voltages[DimSize(voltages,0)-1], "V", cube_3D
cube_3D=Nan
//populate 3D cube from 2D slices
variable v,i,j
for(v=0;v<DimSize(voltages,0);v+=1) //loop across voltage
string temp_str=stringfromlist(v,slices,";")
print temp_str
wave temp_wv=$temp_str
for(i=0;i<DimSize(cube_3D,0);i+=1) //loop across k
for(j=0;j<DimSize(cube_3D,1);j+=1) // loop across energy
cube_3D[i][j][v]=temp_wv[q][p]
endfor
endfor
endfor
//#########INTERPOLATION SECTION#########-> nb the code works totally fine up to this point
//interpolate to fill in the missing data slices
//for(i=0;i<DimSize(cube_3D,0);i+=1) //loop across k
for(i=10;i<20;i+=1) //loop across k
//extract a temp 2D slice in the v,E plane for a given k index
Make/O/N=(DimSize(cube_3D,2),DimSize(cube_3D,1)) temp_2D
SetScale/P x, DimOffset(cube_3D,2), DimDelta(cube_3D,2), "Voltage", temp_2D
SetScale/P y, DimOffset(cube_3D,1), DimDelta(cube_3D,1), "eV", temp_2D
temp_2D=cube_3D[i][q][p]
//NewImage/K=0/F temp_2D //display for debugging purposes
//creat a temp wave to write interpolate values to
duplicate/O temp_2D, interp_temp_2D
//populate iterpolation wave
for(v=0;v<DimSize(voltages,0);v+=1) // loop across voltage
for(j=0;j<DimSize(cube_3D,1);j+=1) // loop across energy
Variable voltage_value, energy_value
voltage_value=DimOffset(cube_3D, 2)+DimDelta(cube_3D,2)*v //work out voltage (ie z) value, not index
energy_value=DimOffset(cube_3D, 1)+DimDelta(cube_3D,1)*j // work out energy (ie y) value, not index
interp_temp_2D=Interp2D(temp_2D, voltage_value, energy_value) // its this bit that doesn't work
endfor
endfor
//NewImage/K=0/F interp_temp_2D //display for debugging purposes
endfor
Killwaves/Z voltages, k_scale, temp_wv, temp_7, temp_neg6, temp_neg4, temp_neg2, temp_1, temp_3, temp_5
End Function
//make waves to hold the information for the slices and correponding volatges
string slices="SMP1212_028_kx;temp_neg6;SMP1212_042_kx;temp_neg4;SMP1212_041_kx;temp_neg2;SMP1212_039_kx;SMP1212_037_kx;temp_1;SMP1212_036_kx;temp_3;SMP1212_035_kx;temp_5;SMP1212_034_kx;temp_7;SMP1212_033_kx"
Make/O voltages={-7, -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8}
String k_scale_str=stringfromlist(0,slices,";")
wave k_scale = $k_scale_str
//make temp waves to fill in missing voltages with NaN waves
duplicate/O k_scale, temp_7
temp_7=nan
duplicate/O temp_7, temp_neg6, temp_neg4, temp_neg2, temp_1, temp_3, temp_5
//make 3D cube to poulate
Make/O/N=(DimSize(k_scale,1), DimSize(k_scale,0), DimSize(voltages,0)) cube_3D
SetScale/P x, DimOffset(k_scale,1), DimDelta(k_scale,1), "Å^-1", cube_3D
SetScale/P y, DimOffset(k_scale,0), DimDelta(k_scale,0), "eV", cube_3D
SetScale/I z, voltages[0], voltages[DimSize(voltages,0)-1], "V", cube_3D
cube_3D=Nan
//populate 3D cube from 2D slices
variable v,i,j
for(v=0;v<DimSize(voltages,0);v+=1) //loop across voltage
string temp_str=stringfromlist(v,slices,";")
print temp_str
wave temp_wv=$temp_str
for(i=0;i<DimSize(cube_3D,0);i+=1) //loop across k
for(j=0;j<DimSize(cube_3D,1);j+=1) // loop across energy
cube_3D[i][j][v]=temp_wv[q][p]
endfor
endfor
endfor
//#########INTERPOLATION SECTION#########-> nb the code works totally fine up to this point
//interpolate to fill in the missing data slices
//for(i=0;i<DimSize(cube_3D,0);i+=1) //loop across k
for(i=10;i<20;i+=1) //loop across k
//extract a temp 2D slice in the v,E plane for a given k index
Make/O/N=(DimSize(cube_3D,2),DimSize(cube_3D,1)) temp_2D
SetScale/P x, DimOffset(cube_3D,2), DimDelta(cube_3D,2), "Voltage", temp_2D
SetScale/P y, DimOffset(cube_3D,1), DimDelta(cube_3D,1), "eV", temp_2D
temp_2D=cube_3D[i][q][p]
//NewImage/K=0/F temp_2D //display for debugging purposes
//creat a temp wave to write interpolate values to
duplicate/O temp_2D, interp_temp_2D
//populate iterpolation wave
for(v=0;v<DimSize(voltages,0);v+=1) // loop across voltage
for(j=0;j<DimSize(cube_3D,1);j+=1) // loop across energy
Variable voltage_value, energy_value
voltage_value=DimOffset(cube_3D, 2)+DimDelta(cube_3D,2)*v //work out voltage (ie z) value, not index
energy_value=DimOffset(cube_3D, 1)+DimDelta(cube_3D,1)*j // work out energy (ie y) value, not index
interp_temp_2D=Interp2D(temp_2D, voltage_value, energy_value) // its this bit that doesn't work
endfor
endfor
//NewImage/K=0/F interp_temp_2D //display for debugging purposes
endfor
Killwaves/Z voltages, k_scale, temp_wv, temp_7, temp_neg6, temp_neg4, temp_neg2, temp_1, temp_3, temp_5
End Function
If you do determine that you need to interpolate, it seems to me that your example data could be interpolated simply by averaging the two adjacent Voltage columns. Interp2d() is not going to be helpful here because you have NaNs in the immediate neighborhood of each interpolated value. My advice is to use a simple MatrixOP command for every NaN column (say colC) so that:
MatrixOP/O interpolatedCol=0.5*(col(temp_2d,colC-1)+col(temp_2d,colC+1))
If your voltage values are in rows, not columns then replace the col() above with row().
A.G.
WaveMetrics, Inc.
October 25, 2017 at 08:26 am - Permalink