The algorithm for making a surface from a set of points in 3D space is to interpolate a grid from that set of points. For example, if the X range of your set of points is -130 to -140 and the Y range of your points is +52 to +54, then one might specify that a grid of 400 Z values might be interpolated for X coordinates that increment by 0.5 units starting at -130 and Y coordinates that increment by 0.1 units starting at +52. This will give you the shape of your surface.
Whether you can color this surface meaningfully now depends on whether you also can meaningfully interpolate the fourth column color values in your grid. Looking at the small amount of example data in your original message, it is hard for me to see what the color correlates to. If it correlates to X/Y-location, then the task is simple; you can interpolate this 4th dimension to X and Y as above, and easily tell IDL to color each vertex of the grid described above with the value calculated by this second interpolation. If the color correlates just to the Z value, then you can use the output of the interpolation in the first paragraph instruction above to determine the color. (In either case, you may need to figure out how to map floating point results for your color to 8-bit integer color values like 40, 50, 60, etc.).
In any case, that is the discussion of the general algorithm. The specifics of the IDL language are that you probably need to do this gridding work outside of ISURFACE before you call iSurface. Thus, you can import the data as Abduwasit Ghulam suggested:
data = READ_ASCII(datafile)
X = data.field1[0,*]
Y = data.field1[1,*]
Z = data.field1[2,*]
D = data.field1[3,*]
but before you pass X, Y and Z to ISURFACE, you would need to do the gridding work outside of ISURFACE (which normally does this for you). The reason is that you probably need the results of the X, Y and Z gridding before you can calculate a 'griddedD' that ISURFACE can use. The easiest IDL command for gridding is GRIDDATA. Thus, here is the GRIDDATA call for the scenario I described in paragraph 1 above:
gridX = -130.0 + (findgen(20) * 0.5)
gridY = 52.0 + (findgen(20) * 0.1)
griddedZ = griddata(x, y, z, /GRID, XOUT=gridX, YOUT=gridY)
; Here is the algorithm if D correlates with X/Y location. The BYTE call will
; convert (by truncating) the floating point result into a color index value.
griddedD = byte(griddata(x, y, z, /GRID, XOUT=gridX, YOUT=gridY))
; Now you need to make your color table
myCT = bytarr(256,3)
; This is just an example of setting one index to red
myCT[40, *] = [255, 0, 0] ; ... etc.
; You would need a color assigned for every value that can appear in the
; 'griddedD' result.
isurface, griddedZ, gridX, gridY, VERT_COLORS=griddedD, RGB_TABLE=myCT
James Jones
|