// $BACKGROUND$
//***Large Dynamic Range Diffraction Pattern***\\
//** Richard Beanland r.beanland@warwick.ac.uk **\\
//* Coding started March 2011 *\\

// Version 1.5, 2 November 2011
// Changes from 1.2: added variable exposure factor
// Changes from 1.3: made script more concise; changed some vasriable names;
// 1.4 - use function Stitch to merge images
// 1.5 - corrected bugs in Stitch

//Define Functions//
//FUNCTION remove outliers
//Remove outliers by comparing with an identical image with median filer applied
void RemoveOutliers(image IMG, number thr)
{
number t,l,b,r
image medImg
GetSelection(IMG,t,l,b,r)
medImg := RealImage("medImg",4,r,b)
medImg=medianFilter(IMG,3,1)
// replace ONLY those pixels in IMG which are >thr in comparison with median image
IMG=tert( (abs(IMG-medImg)>thr), medImg, IMG)
//Delete median image
DeleteImage(medImg)
}//end of RemoveOutliers

//FUNCTION Stitch
image Stitch(image Img1, image Img2, image mask, number sat)
{
number scale
//image mask1
//mask1:=tert( (mask==1),1,0)
//result("Max Img1 = "+max(Img1)+"\n")
//result("Max Img2 = "+max(Img2)+"\n")
//use to find scaling factor to stitch images together
scale=sum((mask==1)*Img1)/sum((mask==1)*Img2)
result("Scaling factor = "+sum((mask==1)*Img1)+"/"+sum((mask==1)*Img2)+" = "+scale+"\n")
Img1=tert( (mask==2),Img2*scale,Img1 )
return Img1
}//End of function Stitch

//FUNCTION Acquire SADP
image AcquireSADP(image Img1)
{
//Define variables
number binning,thresh,flag,mult
number expmax,expo,i,stringlength,t,l,b,r,sat,satHi
string name,pathname,nameImg
image Img2,mask

// USER INPUT
expmax=60
binning=1
thresh=50
sat=10000
mult=2
nameImg = "Diffraction Pattern"
name=" "
// prompt input for File name and location
if (!GetString("File name?",nameImg,nameImg))
exit(0)
if (!SaveAsDialog("","Select the save location - then click Save", pathname)) exit(0)
stringlength=len(pathname)
pathname=left(pathname,stringlength-42)
// prompt input for max exposure time for each image
if (!GetNumber("Maximum exposure time(sec)?",expmax,expmax))
exit(0)
// prompt input for factor difference in exposure time
if (!GetNumber("Factor for exposure time(sec)?",mult,mult))
exit(0)
// prompt input for outlier threshold value
//if (!GetNumber("Threshold for outliers?",thresh,thresh))
//exit(0)
// prompt input for saturation threshold
//if (!GetNumber("Saturation value?",sat,sat))
//exit(0)
// prompt input for binning of images
while (2<3)
{
if (!GetNumber("Image Binning - 1, 2, 4, 8x ?",binning,binning)) exit(0)
if (binning==1 || binning==2 || binning ==4 || binning==8) break
}

// Image acquisition starts here
t = 0
l = 0
b = 2672/binning //values are for Orius camera
r = 2688/binning
expo=expmax
flag=0
i=1
//Img1 is the final high dynamic range image
result("Acquiring first image, exposure time "+expo+" sec\n")
Img1 := SSCGainNormalizedBinnedAcquire(expo, binning, t, l, b, r)
Img1.RemoveOutliers(thresh)
Img1.ImageSetName(nameImg)
Img1.DisplayAt(150, 40)// Display the image on the screen
if ( max(Img1)<sat )
{
flag=1
}
result("Max Img1 = "+max(Img1)+"\n")
result("saturation = "+sat+"\n")
mask:=tert( ((Img1>0.9*sat)*(Img1<sat)), 1, 0 )
mask=mask+tert( (Img1>sat), 2, 0 )
//mask.DisplayAt(200,50)
mask.ImageSetName("Mask")
expo=expo/mult//reduce exposure time
while (flag==0)//loop until acquired image is not saturated
{
i=i+1
result("Acquiring image "+i+", exposure time "+expo+" sec\n")
Img2 = SSCGainNormalizedBinnedAcquire(expo, binning, t, l, b, r)
Img2.RemoveOutliers(thresh)
Img1.Stitch(Img2,mask,sat)
Img1.UpdateImage()
if ( max(Img2)<sat )
{
flag=1
}
else
{
mask=tert( ((Img2>0.9*sat)*(Img2<sat)), 1, 0)
if (sum(mask)==0)
{
throw ("Mask empty!")
exit(0)
}
mask=mask+tert( (Img2>sat), 2, 0 )
//mask.UpdateImage()
result("Image is still saturated\n")
}
expo=expo/mult//reduce exposure time\
if (expo<0.05)
{
//flag=1
}
}
Img2.DeleteImage()
result("\nImage Max,Min = "+max(Img1)+", "+min(Img1)+"\n")
nameImg=pathname+nameImg
SaveasGatan(Img1,nameImg)//Save the image to a folder
return Img1
}//end of acquireSADP

//************END OF FUNCTIONS - MAIN PROGRAM FOLLOWS ******************
//It's not very long...
image dynImg=AcquireSADP(dynImg)
result("Acquisition complete\n")
result("You may now lower the screen\n\n")