Index: /trunk/stac/scripts/mkfakeimages.pl
===================================================================
--- /trunk/stac/scripts/mkfakeimages.pl	(revision 742)
+++ /trunk/stac/scripts/mkfakeimages.pl	(revision 742)
@@ -0,0 +1,145 @@
+#!/usr/local/bin/perl
+#
+# testcombine.pl
+#
+# Program to create some dithered/rotated frames and combine them
+# to see how the combination works.
+#
+#
+# Paul A. Price
+#
+# CVS: $Id: mkfakeimages.pl,v 1.1 2004-05-20 00:06:21 price Exp $
+#
+# History:
+# Added orderly star placement.  Fixed so that the star positions are
+# accurate once mapped.  --- PAP, 20/10/03.
+# Added FITS headers to the created images so that the combination
+# program uses the correct statistics.  --- PAP 22/10/03.
+#
+
+
+# Parameters.
+$DEG_IN_RAD = 57.2957795130823;   # Degrees per radian.
+
+$GAIN = 1;			# Number of electrons per ADU
+$READNOISE = 10;		# Detector read noise (electrons);
+$NUMSTARS = 500;		# Number of stars on images.
+$INSIZE = 1024;		    # Size of input images (pixels on a side).
+$OUTSIZE = 1024;	    # Size of output image (pixels on a side).
+$SEEING = 0; # Seeing in input images (input pixels).  Zero means use array.
+$NUMCRS = 0;			# Number of CRs on input images.
+$STARFILE = "stars";	      # File containing list of stars (x,y,f).
+$IMAGE = "test";		# Base file name for images and stuff.
+$BACKGROUND = 10000;		# Sky background for input images
+$OUTSCALE = sqrt(1);  # Pixel scale of output image compared to input.
+$STARANGLE = 1.5; # Angle at which a regular star field is tilted (degrees)
+$STARFLUXMIN = 500;		# Minimum flux of stars (ADU)
+$STARFLUXMAX = 20000;		# Maximum flux of stars (ADU)
+$MAXOFFSET = 5;	 # Maximum offset for each image (output image scale).
+@ROTATIONS=(0,3);		# Rotations of the images
+@SEEING = (2.5, 4.5);	       # Seeing in input images (input pixels)
+
+
+$CREATE = "./fakeimage -d -x $INSIZE -y $INSIZE -g $GAIN -r $READNOISE -b $BACKGROUND";   # Command to create an input image
+$CREATELIST = "-m ";   # Switch to take stars from a list.
+$CREATESEEING = "-f ";   # Switch for seeing
+
+# Parse command-line arguments
+for ($i=0;$i<=$#ARGV;$i++)
+{
+    if ($ARGV[$i] eq "-numstars") { $NUMSTARS = $ARGV[++$i]; }
+    elsif ($ARGV[$i] eq "-insize") { $INSIZE = $ARGV[++$i]; }
+    elsif ($ARGV[$i] eq "-outsize") { $OUTSIZE = $ARGV[++$i]; }
+    elsif ($ARGV[$i] eq "-seeing") { $SEEING = $ARGV[++$i]; }
+    elsif ($ARGV[$i] eq "-numcrs") { $NUMCRS = $ARGV[++$i]; }
+    elsif ($ARGV[$i] eq "-switches") { $COMBINESWITCHES = $ARGV[++$i]; }
+    else
+    {
+	die "Unknown command-line switch: $ARGV[$i]\n";
+    }
+}
+
+
+
+# Generate list of stars.
+@x=();   # List of star x positions
+@y=();   # List of star y positions
+@f=();   # List of star fluxes
+open STARS, "> $IMAGE.stars";
+$num=0;
+for ($i=0;$i<sqrt($NUMSTARS);$i++)
+{
+    # Random distribution
+#    push @x, rand($OUTSIZE);
+#    push @y, rand($OUTSIZE);
+#    push @f, rand(32*1024-$BACKGROUND);
+
+    # Ordered distribution
+    for ($j=0;$j<sqrt($NUMSTARS);$j++)
+    {
+	push @x, ($i+1)*$OUTSIZE/(sqrt($NUMSTARS)+1)*cos($STARANGLE/$DEG_IN_RAD) + ($j+1)*$OUTSIZE/(sqrt($NUMSTARS)+1)*sin($STARANGLE/$DEG_IN_RAD);
+	push @y, - ($i+1)*$OUTSIZE/(sqrt($NUMSTARS)+1)*sin($STARANGLE/$DEG_IN_RAD) + ($j+1)*$OUTSIZE/(sqrt($NUMSTARS)+1)*cos($STARANGLE/$DEG_IN_RAD);
+	if ($STARFLUXMIN == $STARFLUXMAX) { push @f, $STARFLUXMAX; }
+	else { push @f, rand($STARFLUXMAX-$STARFLUXMIN) + $STARFLUXMIN; }
+
+	print STARS "$x[$num] $y[$num] $f[$num]\n";
+	$num++;
+    }
+}
+
+
+# Iterate over the input images
+for ($i=0;$i<=$#ROTATIONS;$i++)
+{
+    # Get transformation.
+    $B00 = cos($ROTATIONS[$i]/$DEG_IN_RAD) / $OUTSCALE;
+    $B01 = sin($ROTATIONS[$i]/$DEG_IN_RAD) / $OUTSCALE;
+    $B10 = - $B01;
+    $B11 = $B00;
+    $A0 = $INSIZE/2;
+    $A1 = $INSIZE/2;
+    if ($MAXOFFSET == 0)
+    {
+	$offsetx = 0;
+	$offsety = 0;
+    }
+    else
+    {
+	$offsetx = rand(2*$MAXOFFSET) - $MAXOFFSET;
+	$offsety = rand(2*$MAXOFFSET) - $MAXOFFSET;
+    }
+
+    # Output star positions to file
+    open STARS, "> ${IMAGE}_$i.stars";
+    for ($j=0;$j<$num;$j++)
+    {
+	$x = $B00*($x[$j] - 1.0 - $OUTSIZE/2 - $offsetx) + $B10*($y[$j] - 1.0 - $OUTSIZE/2 - $offsety) + $A0;
+	$y = $B01*($x[$j] - 1.0 - $OUTSIZE/2 - $offsetx) + $B11*($y[$j] - 1.0 - $OUTSIZE/2 - $offsety) + $A1;
+
+
+	print STARS "$x $y $f[$j]\n";
+    }
+    close STARS;
+
+
+    # Create image
+    if ($SEEING == 0) { $seeing = $SEEING[$i]; }
+    else { $seeing = $SEEING; }
+    print "$CREATE $CREATESWITCHES $CREATELIST ${IMAGE}_$i.stars $CREATESEEING $seeing ${IMAGE}_$i.fits\n";
+    system "$CREATE $CREATESWITCHES $CREATELIST ${IMAGE}_$i.stars $CREATESEEING $seeing ${IMAGE}_$i.fits";
+
+    # Create map
+    $B00 = $OUTSCALE*cos($ROTATIONS[$i]/$DEG_IN_RAD);
+    $B01 = - $OUTSCALE*sin($ROTATIONS[$i]/$DEG_IN_RAD);
+    $B10 = - $B01;
+    $B11 = $B00;
+
+    $A0=$OUTSIZE/2 - ($B00 + $B10)*$INSIZE/2 + $offsetx;
+    $A1=$OUTSIZE/2 - ($B01 + $B11)*$INSIZE/2 + $offsety;
+
+    # Output map
+    open MAP, "> ${IMAGE}_$i.fits.map";
+    print MAP "$A0 $A1\n$B00 $B01\n$B10 $B11\n";
+    close MAP;
+}
+