RGB color stack

// $BACKGROUND$

image dataset := GetFrontImage()
number roff,goff,boff
number wavelength 

roff = 1
goff = .5
boff = .2

Number Nx = ImageGetDimensionSize(dataset, 0)
Number Ny = ImageGetDimensionSize(dataset, 1)
Number Nz = ImageGetDimensionSize(dataset, 2)
//Result("\n number of slices: " + Nz )

image NewData := RGBimage("data",4, Nx, Ny, Nz)
image rawimage := Realimage("grayframe",4, Nx, Ny)
image newimage := RGBimage("singleframe",4, Nx, Ny, 1)
number Scale = dataset.ImageGetDimensionScale(2)
number Origin = dataset.ImageGetDimensionOrigin(2)
string Unit = dataset.ImageGetDimensionUnitString(2)
number Scaledeg = dataset.ImageGetDimensionScale(1)
number Origindeg = dataset.ImageGetDimensionOrigin(1)
string Unitdeg = dataset.ImageGetDimensionUnitString(1)
number sliceTowavelength(number plane, image data)
{

	wavelength = Scale*plane + Origin
	return wavelength
	
}


number gamma = 0.8 // display adjustment for monitor type
number attenuation = 1

 
number j
for (j=0; j<Nz; j++)

{




	wavelength = sliceToWavelength(j,dataset)

	If (wavelength < 360)
	{
	roff = 1
	goff = 1.0 
	boff = 1
	
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB(roff*rawimage,goff*rawimage,boff*rawimage)
	} 

	If (wavelength >= 360 && wavelength < 380)
	{
	roff = 1
	goff = ((380-wavelength)/(380-360))**gamma
	boff = 1.0
	
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB(roff*rawimage,goff*rawimage,boff*rawimage)
	} 


	// if you can get wavelength:
	If (wavelength >= 380 && wavelength <= 440)	
	{
	attenuation = 0.3 + 0.7 * (wavelength - 380) / (440 - 380) 
	roff = ((-(wavelength - 440) / (440 - 380)) * attenuation)**gamma + 1 - (1.0 * attenuation)**gamma 
	goff = 0.0 
	boff = 1.0
	//boff = (1.0 * attenuation)**gamma 
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB((roff)*rawimage,goff*rawimage,(boff)*rawimage)
	} 
	
	If (wavelength > 440 && wavelength <= 490)
	{
	//attenuation = 0.3 + 0.7 * (wavelength - 380) / (440 - 380) 
	roff = 0
	goff = ((wavelength-440)/(490-440))**gamma 
	boff = 1 
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB(roff*rawimage,goff*rawimage,boff*rawimage)
	} 
	
	If (wavelength > 490 && wavelength <= 510)
	{
	//attenuation = 0.3 + 0.7 * (wavelength - 380) / (440 - 380) 
	roff = 0
	goff = 1.0 
	boff = ((510-wavelength)/(510-490))**gamma
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB(roff*rawimage,goff*rawimage,boff*rawimage)
	}
	
	If (wavelength > 510 && wavelength <= 580)
	{
	//attenuation = 0.3 + 0.7 * (wavelength - 380) / (440 - 380) 
	roff = (((wavelength - 510) / (580 - 510)))**gamma 
	goff = 1
	boff = 0
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB(roff*rawimage,goff*rawimage,boff*rawimage)
	}
	
	If (wavelength > 580 && wavelength <= 645)
	{
	//attenuation = 0.3 + 0.7 * (wavelength - 380) / (440 - 380) 
	roff = 1
	goff = (-(wavelength-645)/(645-580))**gamma 
	boff = 0
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB(roff*rawimage,goff*rawimage,boff*rawimage)
	}

	If (wavelength >= 645 && wavelength <= 750)
	{
	attenuation = 0.3 + 0.7 * (750 - wavelength) / (750 - 645) 
	roff = 1 
	goff = 1*(((wavelength - 645) / (750 - 645)))**gamma
	boff = 1*(((wavelength - 645) / (750 - 645)))**gamma
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB(roff*rawimage,goff*rawimage,boff*rawimage)
	}

	If (wavelength >= 750)
	{
	//attenuation = 0.3 + 0.7 * (750 - wavelength) / (750 - 645) 
	roff = 1
	goff = 1 
	boff = 1
	rawimage = Slice2(dataset, 0, 0, j, 0, Nx, 1, 1, Ny, 1)
	rawimage = 255*rawimage/max(rawimage) //4.76237e7
	newimage = RGB(roff*rawimage,goff*rawimage,boff*rawimage)
	}
	

	
	Slice2(NewData, 0, 0, j, 0, Nx, 1, 1, Ny, 1) = newimage
}
Result("\n Scale: " + Scale + " " + Unit + " /frame " + "; & Origin: " + Origin + " " + Unit)
ImageSetDimensionScale(NewData,2,Scale)
ImageSetDimensionOrigin(NewData,2,Origin)
ImageSetDimensionUnitString(NewData,2,Unit)

ImageSetDimensionScale(NewData,0,Scaledeg)
ImageSetDimensionOrigin(NewData,0,Origindeg)
ImageSetDimensionUnitString(NewData,0,Unitdeg)

ImageSetDimensionScale(NewData,1,Scaledeg)
ImageSetDimensionOrigin(NewData,1,Origindeg)
ImageSetDimensionUnitString(NewData,1,Unitdeg)
NewData.showImage()