/Virgo/FeldPNe/PNe8stattest/PNeNstat.sm
OR
/Virgo/FeldPNe/vat/PNe8stattest/PNeNstat.sm
PNeNstat
define TeX_strings 1
#program to generate graph of # pixels counted twice
as a function of number of boxes
#will need to generate .wfound.dat .Rfound.dat
.0found.dat files for
# rang : 10 to 60 by 5s
# N boxes : 25 to 200 by 25s
# imag : 3,3m,4,4m,7,7m,Sub,Subm
#naming method will be:
#8.50.Im3.fits.R10.wfound.dat
# where
# 8 is version number
# 50 is N_boxes
# Im3.fits is imag = 3.fits
# R10 is rang = 10
# wfound means boxes each have at
least 1 PN
#establish running parameters
set ranges = 10, 60, 5
set images = { 3.fits 3m.fits 4.fits 4m.fits 7.fits
7m.fits Sub.fits Subm.fits }
set Nbs = 25, 200, 25
#for testing
#set images = { 4.fits }
#set ranges = { 10 15 20 }
#set Nbs = { 25 50 75 100 }
#add to normal colours
add_ctype orange 255 153 0
add_ctype purple 153 0 255
add_ctype brown 102 51 0
add_ctype grey 204 204 204
add_ctype perriwinkle 153 153 255
add_ctype olive 102 102 0
add_ctype blood 153 0 0
set Colours = { black red green blue cyan magenta
orange purple brown perriwinkle olive blood }
data v.v
read { vnum 1 }
define v $(vnum[0])
data PNe.param
read { SZ 1 SO 2 SR 3 }
define SpaceZero $(SZ[0])
define SpaceOne $(SO[0])
#define SpaceTwo $(ST[0])
define SpaceRand $(SR[0])
foreach imag images {
define file_type fits
image $imag
PNew
foreach rang ranges {
define rangeRD
$($($rang)*ABS($CD1_1))
#call max
finding program to find max # pixels, write to file
# ONLY
if necessary!!!
#PNemax
#temp! delete
PNemaxread line here later
PNemaxread
#if maxes have
been computed, can use this one, but don't have to
#PNemaxread
foreach Nb Nbs
{
#run .dat file generator
# this is basically a (slightly
modified) copy of PNe.sm
#ONCE this has been run for everything, comment out
and only re-do graphing/pixel-finding for speed
#PNeNdat
}
}
#run doubled pixel finder/grapher
#if this has already run, only
run the program PNeG-->
#PNeD
PNeG
#PNeG can only be run if 8.N.etc
has been done for selected images/ranges
}
PNeNdat
#generate all .dat files for selected image, range,
and N_boxes
set DIMEN(wPNe) = 0
set DIMEN(woPNe) = 0
set DIMEN(rand) = 0
define reals $Nb
echo running image = $imag range =
$rang N_boxes = $Nb
#run boxes with PNe
PNe1
#run boxes without PNe
PNe0
#run random boxes
PNer
PNew
#this does the real work and is called by PNe, which
is just loops
#define as fits and read in important header info
#x-axis length
define NAXIS1 image
#y-axis length
define NAXIS2 image
#0 point in RA/DEC
define CRVAL1 image
define CRVAL2 image
#0 point in image x/y
define CRPIX1 image
define CRPIX2 image
#pixel scale in dRA/dx and dDEC/dy
define CD1_1 image
define CD2_2 image
#prep for graphing inline
set bins = -5, 100, 0.5
#read in PNe datafile
data PNe.dat
read { field 1 num 2 RAh 3 RAm 4 RAs 5 DECd 6 DECam
7 DECas 8 m5007 9 notes 10.s }
#set RA/DEC both in degrees and range, too.
set DEC = DECd + (DECam/60) + (DECas/3600)
set RA = 15 * (RAh + (RAm/60) + (RAs/3600))
#^these conversions are checked and good. duh.^#
PNe1
#random boxes that contain one (or more) PN
echo starting boxes WITH PNe...
#file to save box coordinates in and compare to later
define wfound "$!v.$!Nb.Im$!imag.R$!rang.wfound.dat"
set dimen(saveRAw) = $reals
set dimen(saveDECw) = $reals
#loop for as many boxes as we want
define makesure 0
define num 0
while {$num < $reals} {
#generate random (x,y) coordinates
define rx $(($(RANDOM(1))*($NAXIS1-(2*$rang)))+$rang)
define ry $(($(RANDOM(1))*($NAXIS2-(2*$rang)))+$rang)
#convert to RA/DEC
define rDEC (($ry - $CRPIX2 + 1)*$(ABS($CD2_2)) + $CRVAL2)
define rRA $($(-($rx - $CRPIX1 + 1 ))*$(ABS($CD1_1)/cosd($rDEC)) +
$CRVAL1)
#looking for PNe in a box
centered around the point that is 2*rang x 2*range
set delRA=$rRA-RA
set delDEC=$rDEC-DEC
set counter=( (abs(delRA)<=$rangeRD) &&
(abs(delDEC)<=$rangeRD) )? 1 : 0
if
($(sum(counter)) > 0) {
#found at
least one PN in this box
#check for
other boxes too close, if SpaceOne is one
set delRA =
$rRA - saveRAw
set delDEC = $rDEC - saveDECw
#check
for SpaceOne
if ($SpaceOne == 1) {
set counterp = ( (abs(delRA)<=$($rangeRD)) &&
(abs(delDEC)<=$rangeRD) )? 1 : 0
} ELSE {
set counterp = {0}
}
if
($(sum(counterp)) == 0) {
#now we have a box that contains 0 and is far enough from other boxes
#save coords and save pixels
set saveRAw[$num] = $rRA
set saveDECw[$num] = $rDEC
define num $($num + 1)
# don't need to add up
pixels, only concerned with boxes, for counting how many times each
pixel gets counted
#
set dimen(statra)=$($(2*$rang+1)*$(2*$rang+1))
#
define idum 0
#
do x = -($rang), $rang {
#
do y = -($rang), $rang {
#
set
statra[$idum]=$(image[$rx+$x,$ry+$y])
#
#handles masked pixels
#
set statrafinal=statra if (statra>-900)
#
define idum $($idum+1)
#
}
#
}
#
#
#concat this image onto the growing big file
#
set wPNe = wPNe concat statrafinal
#end of
checking that other boxes aren't too close
} ELSE {
define makesure $($makesure + 1)
}
#end of is 1+ PN found
} ELSE {
define
makesure $($makesure + 1)
}
#check makesure - don't let it
look more than 10000 times #_ boxes
if ($makesure >
$($reals*10000)) {
echo There is trouble - Couldn't get enough boxes in $($reals*10000)
trials... you suck
write "$!v.$!Nb.Im$!imag.R$!rang.problems.error.wfound.z" too many
loops in wfound - 1 PNe
define numz
$num
define num $reals
}
#end while still doing more boxes, while $num <
$reals, loop
}
#write out data file
!rm $wfound
do w = 0,
$(dimen(saveRAw)-1) { write $wfound "$!(saveRAw[$!w])
$!(saveDECw[$!w]) $!rang" }
echo 1/3 - did 1+ PNe/box
WITH optional spacing
PNe0
#checking for boxes WITHOUT PNe in them. # boxes =
$reals
echo starting boxes WITHOUT PNe
#file to save box coordinates in and compare to later
define zero "$!v.$!Nb.Im$!imag.R$!rang.0found.dat"
set dimen(saveRAz) = $reals
set dimen(saveDECz) = $reals
define makesure 0
define num 0
while {$num < $reals} {
#pick random (x,y)
define rx
$(($(RANDOM(1))*($NAXIS1-(2*$rang)))+$rang)
define ry
$(($(RANDOM(1))*($NAXIS2-(2*$rang)))+$rang)
#convert to RA/DEC
define rDEC (($ry - $CRPIX2
+ 1)*$(ABS($CD2_2)) + $CRVAL2)
define rRA $($(-($rx - $CRPIX1 +
1 ))*$(ABS($CD1_1)/cosd($rDEC)) + $CRVAL1)
#looking for PNe in a box
centered around the point that is 2*rang x 2*range
set delRA=$rRA-RA
set delDEC=$rDEC-DEC
set counter=(
(abs(delRA)<=$rangeRD) && (abs(delDEC)<=$rangeRD) )? 1 : 0
if
($(sum(counter)) == 0) {
#make sure no
other boxes close to here
set delRA =
$rRA - saveRAz
set delDEC =
$rDEC - saveDECz
#check for
SpaceZero
if ($SpaceZero
== 1) {
set counterp = ( (abs(delRA)<=$($rangeRD))
&& (abs(delDEC)<=$rangeRD) )? 1 : 0
} ELSE {
set counterp = {0}
}
if
($(sum(counterp)) == 0) {
#now we have a box that contains 0 and is far enough
from other boxes
set saveRAz[$num] = $rRA
set saveDECz[$num] = $rDEC
define num $($num + 1)
#
#collecting from a box around PNe 2*rang+1 x 2*rang+1
#
#
set dimen(statra)=$($(2*$rang+1)*$(2*$rang+1))
#
define idum 0
#
do x = -($rang), $rang {
#
do y = -($rang), $rang {
#
set
statra[$idum]=$(image[$rx+$x,$ry+$y])
#
#handles
masked pixels
#
set
statrafinal=statra if (statra>-900)
#
define idum
$($idum+1)
#
}
#
}
#
#
#concat this image onto the growing big file
#
set woPNe = woPNe concat statrafinal
#end of making
sure no other boxes too close
} ELSE {
define makesure $($makesure + 1)
}
#end of if no PNe found
} ELSE {
define
makesure $($makesure + 1)
}
#check to make sure it doesn't
loop too long
if ($makesure >
$($reals*1000)) {
echo There is
trouble - Couldn't get enough boxes in $($reals*10000) trials... you
suck
write "$!v.$!Nb.Im$!imag.R$!rang.problems.error.0found.z" too many
loops in 0found - 0 PNe
define numz
$num
define num
$reals
}
#end while $num < $reals loop for overall 0 PNe
loop
}
#write out data file
!rm $zero
do w = 0, $(dimen(saveRAz)-1) { write $zero
"$!(saveRAz[$!w]) $!(saveDECz[$!w]) $!rang" }
echo 2/3 - did 0 PNe/box WITH optional
spacing
PNer
#random boxes
echo starting random boxes...
#initialise storage vars/vects
set dimen(saveRAr) = $reals
set dimen(saveDECr) = $reals
#file to save box coordinates in and compare to later
define r "$!v.$!Nb.Im$!imag.R$!rang.Rfound.dat"
define numr 0
define makesurer 0
while {$numr < $reals} {
#generate random coordinates as
before
define rx
$(($(RANDOM(1))*($NAXIS1-(2*$rang)))+$rang)
define ry
$(($(RANDOM(1))*($NAXIS2-(2*$rang)))+$rang)
#convert to RA/DEC
define rDEC (($ry - $CRPIX2
+ 1)*$(ABS($CD2_2)) + $CRVAL2)
define rRA $($(-($rx - $CRPIX1 +
1 ))*$(ABS($CD1_1)/cosd($rDEC)) + $CRVAL1)
#make sure no other boxes close
to here
set delRA = $rRA - saveRAr
set delDEC = $rDEC - saveDECr
#check for SpaceRand
if ($SpaceRand == 1) {
set counterp =
( (abs(delRA)<=$($rangeRD)) && (abs(delDEC)<=$rangeRD) )?
1 : 0
} ELSE {
set counterp =
{0}
}
if
($(sum(counterp)) == 0) {
set
saveRAr[$numr] = $rRA
set
saveDECr[$numr] = $rDEC
define numr
$($numr + 1)
#collecting
from a box around PNe 2*rang+1 x 2*rang+1
# set
dimen(statra)=$($(2*$rang+1)*$(2*$rang+1))
# define idum 0
# do x =
-($rang), $rang {
#
do y = -($rang), $rang {
#
set
statra[$idum]=$(image[$rx+$x,$ry+$y])
#
#handles masked pixels
#
set statrafinal=statra if
(statra>-900)
#
define idum $($idum+1)
#
}
#
# }
#
# #concat this
image onto the growing big file
# set rand =
rand concat statrafinal
#end of if sum == 0 loop
} ELSE {
define
makesurer $($makesurer + 1)
}
if ($makesurer >
$($reals*10000)) {
echo There is
trouble - Couldn't get enough boxes in $($reals*10000) trials... you
suck
write "$!v.$!Nb.Im$!imag.R$!rang.problems.error.Rfound.z" too many
loops in Rfound - you must have things set very wrong
define numrr
$numr
define numr
$reals
}
#end of checking through random boxes
}
#write data file out
!rm $r
do wr = 0, $(dimen(saveRAr)-1) { write $r
"$!(saveRAr[$!wr]) $!(saveDECr[$!wr]) $!rang" }
echo 3/3 - did random boxes
PNePG
#prep graphing device
device x11
# device postencap "crap.Im$!imag.ps"
expand .6
erase
ctype black
define x_gutter 0.5
#prep labels on right side of screen
window 4 1 4 1
lim 0 15 0 15
do clr = 0, $(dimen(ranges)-1) {
relocate 0 $clr
ctype $(Colours[$clr])
label range = $(ranges[$clr]) -
$(Colours[$clr])
}
ctype black
vecminmax Nbs xmin xmax
define xmax $($xmax + $($xmax/10))
define xmin $($xmin - $($xmin/10))
#these are used elsewhere - set here
define ymin -0.1
define ymax 3
set xlims = {1 1}
set ylims = {1 1}
set xlims[0] = $xmin
set xlims[1] = $xmax
set ylims[0] = $ymin
set ylims[1] = $ymax
#set ylims[0] = -0.1
#set ylims[1] = 3.1
#prep graphing sides
window 4 1 1 1
lim xlims ylims
box
lim 0 10 0 10
relocate 1 10.1
label 1+ PNe
ylabel Fraction of pixels counted more than once /
total possible pixels
window 4 1 2 1
lim xlims ylims
box
lim 0 10 0 10
relocate 1 10.1
label Random
relocate 5 10.35
label $imag
xlabel Number of Boxes
window 4 1 3 1
lim xlims ylims
box
lim 0 10 0 10
relocate 1 10.1
label 0 PNe
xlabel Number of Boxes
PNeD
#part where we read back in the data files and check
doubled pixels
#this is run once for each image, so we have all
ranges and Nb values to check
#prep graphs
PNePG
define it -1
#strategy:
#foreach range ( each Nb size, read in positions of
boxes/ranges )
#go through each method and box in the list and add
1 to those pixels' values
#
methods = 1+ PNe found, random, 0 PNe found
#files will look like:
"$!v.$!Nb.Im$!imag.R$!rang.wfound.dat"
#and have: RA DEC rang
#loop to extract every N_box value at every range
foreach rang ranges {
#iterate counter for ranges
define it $($it + 1)
define rangeRD
(($rang)*ABS($CD1_1))
#read in max pixel values for
this range/image selection
PNemaxread
#variables are: npix1 npix0 npixR
#SECTION TO DO
1+ PNe FOUND
set dimen(2count1) = dimen(Nbs)
define it2 -1
foreach Nb Nbs {
#iterate
counter for 2founds
define it2
$($it2 + 1)
data "$!v.$!Nb.Im$!imag.R$!rang.wfound.dat"
read { wRA 1
wDEC 2 crap 3.s }
#setup 'image'
to store # times counted for each pixel
image (0,0)
image ($NAXIS1, $NAXIS2)
#loop over the
entire box catalog and set pixels in 'image' to the number of times
they've each been counted
do bl = 0, $(DIMEN(wRA)-1) {
#convert current box center to x,y in pixels
define y $(INT($CRPIX2 +
($($(wDEC[$bl])-$CRVAL2)/$(ABS($CD2_2)))))
define x $(INT($($CRPIX1 - 1 - $($($(wRA[$bl] -
$CRVAL1)*cosd($(wDEC[$bl])))/$(ABS($CD1_1))))))
#loop over every pixel in box
do xadd = -$rang, $rang {
do yadd = -$rang, $rang {
#add 1 to pixel
set
image[$($x+$xadd),$($y+$yadd)] = $($(image[$($x+$xadd),$($y+$yadd)]) +
1)
}
}
#end loop over
box .dat file
}
echo counting
up 2s
#count up how many pixels in image have value > 1
do x = 0,
$($NAXIS1-1) {
do y = 0, $($NAXIS2-1) {
if ($(image[$x,$y]) > 1) {
set
2count1[$it2] = $($(2count1[$it2]) + 1)
#end if 2
}
#end x loop
}
#end y loop
}
#define sum 0
##new method:
average the number of counts for all pixels
#do x = 0,
$($NAXIS1-1) {
# do y = 0, $($NAXIS2-1) {
# define sum $($sum +
$(image[$x,$y]))
# }
#}
#set
2count1[$it2] = $($sum/$($NAXIS1 * $NAXIS2))
#echo wfounds
Nb is $(Nbs[$it2]) $sum $($sum/$($NAXIS1*$NAXIS2))
#write this image to output file for later use
write image
"$!v.$!Nb.Im$!imag.R$!rang.wfound.N.fits"
}
#we are concerned w/ ratio, so
divide # pixels by total possible
set 2count1r =
2count1/$npix1
#add this range's data to the
graph for this image
ctype
$(Colours[$it])
window 4 1 1 1
lim xlims ylims
connect Nbs 2count1r
#SECTION TO DO
randoms
set dimen(2countR) = dimen(Nbs)
define it2 -1
foreach Nb Nbs {
#iterate
counter for 2Rfounds
define it2
$($it2 + 1)
data
"$!v.$!Nb.Im$!imag.R$!rang.Rfound.dat"
read { rRA 1
rDEC 2 crap 3.s }
#setup 'image'
to store # times counted for each pixel
image (0,0)
image
($NAXIS1, $NAXIS2)
#loop over the
entire box catalog
do bl = 0,
$(DIMEN(rRA)-1) {
#convert current box center to x,y in pixels
define y $(INT($CRPIX2 +
($($(rDEC[$bl])-$CRVAL2)/$(ABS($CD2_2)))))
define x $(INT($($CRPIX1 - 1 - $($($(rRA[$bl] -
$CRVAL1)*cosd($(rDEC[$bl])))/$(ABS($CD1_1))))))
#loop over every pixel in box
do xadd = -$rang, $rang {
do yadd = -$rang, $rang {
#add 1 to pixel
set
image[$($x+$xadd),$($y+$yadd)] = $($(image[$($x+$xadd),$($y+$yadd)]) +
1)
}
}
#end loop over
box .dat file
}
echo counting
up 2s
#count up how
many pixels in image have value > 1
do x = 0,
$($NAXIS1-1) {
do y = 0, $($NAXIS2-1) {
if ($(image[$x,$y]) > 1) {
set
2countR[$it2] = $($(2countR[$it2]) + 1)
#end if 2
}
#end x loop
}
#end y loop
}
#define sum 0
##new method:
average the number of counts for all pixels
#do x = 0,
$($NAXIS1-1) {
# do y = 0, $($NAXIS2-1) {
# define sum $($sum +
$(image[$x,$y]))
# }
#}
#set
2countR[$it2] = $($sum/$($NAXIS1 * $NAXIS2))
#echo Rfounds
Nb is $(Nbs[$it2]) $sum $($sum/$($NAXIS1*$NAXIS2))
#write this
image to output file for later use
write image
"$!v.$!Nb.Im$!imag.R$!rang.rfound.N.fits"
}
#we are concerned w/ ratio, so
divide # pixels by total possible
set 2countRr = 2countR/$npixR
#add this range's data to the
graph for this image
ctype $(Colours[$it])
window 4 1 2 1
lim xlims ylims
connect Nbs 2countRr
#SECTION TO DO 0 PNe FOUND
set dimen(2count0) = dimen(Nbs)
define it2 -1
foreach Nb Nbs {
#iterate
counter for 2found0s
define it2
$($it2 + 1)
data
"$!v.$!Nb.Im$!imag.R$!rang.0found.dat"
read { zRA 1
zDEC 2 crap 3.s }
#setup 'image'
to store # times counted for each pixel
image (0,0)
image
($NAXIS1, $NAXIS2)
#loop over the
entire box catalog
do bl = 0,
$(DIMEN(zRA)-1) {
#convert current box center to x,y in pixels
define y $(INT($CRPIX2 +
($($(zDEC[$bl])-$CRVAL2)/$(ABS($CD2_2)))))
define x $(INT($($CRPIX1 - 1 - $($($(zRA[$bl] -
$CRVAL1)*cosd($(zDEC[$bl])))/$(ABS($CD1_1))))))
#loop over every pixel in box
do xadd = -$rang, $rang {
do yadd = -$rang, $rang {
#add 1 to pixel
set
image[$($x+$xadd),$($y+$yadd)] = $($(image[$($x+$xadd),$($y+$yadd)]) +
1)
}
}
#end loop over
box .dat file
}
echo counting
up 2s
#count up how
many pixels in image have value > 1
do x = 0,
$($NAXIS1-1) {
do y = 0, $($NAXIS2-1) {
if ($(image[$x,$y]) > 1) {
set
2count0[$it2] = $($(2count0[$it2]) + 1)
#end if 2
}
#end x loop
}
#end y loop
}
#define sum 0
##new method:
average the number of counts for all pixels
#do x = 0,
$($NAXIS1-1) {
# do y = 0, $($NAXIS2-1) {
# define sum $($sum +
$(image[$x,$y]))
# }
#}
#set
2count0[$it2] = $($sum/$($NAXIS1 * $NAXIS2))
#echo 0founds
Nb is $(Nbs[$it2]) $sum $($sum/$($NAXIS1*$NAXIS2))
#write this
image to output file for later use
write image
"$!v.$!Nb.Im$!imag.R$!rang.0found.N.fits"
}
#we are concerned w/ ratio, so
divide # pixels by total possible
set 2count0r = 2count0/$npix0
#add this range's data to the
graph for this image
ctype $(Colours[$it])
window 4 1 3 1
lim xlims ylims
connect Nbs 2count0r
#write this
data out, too, for my sanity
do wtt = 0, $(dimen(2count1)-1) {
write
"$!v.N.Im$!imag.R$!rang.found.N.dat" $(Nbs[$wtt]) $(2count1r[$wtt])
$(2countRr[$wtt]) $(2count0r[$wtt])
}
#end of ranges loop
}
ctype $(Colours[$it])
#use $co as colour name in titles/labels
PNeG
#just re-read in data and graph it, no more
calculations
#prep graphs
PNePG
define it -1
#reconstructs graphs from data already calculated
completely
foreach rang ranges {
#iterate counter for ranges
define it $($it + 1)
data
"$!v.N.Im$!imag.R$!rang.found.N.dat"
read { Nbs 1 2count1 2 2countR 3
2count0 4 }
#for all graphing
ctype $(Colours[$it])
#for graphing 1+ found
window 4 1 1 1
lim xlims ylims
connect Nbs 2count1
#graphing random
window 4 1 2 1
lim xlims ylims
connect Nbs 2countR
#graphing 0 found
window 4 1 3 1
lim xlims ylims
connect Nbs 2count0
}
PNemaxread
#code to read in max values
data "$!v.max.Im$!imag.R$!rang.max"
read { npix1v 1 npix0v 2 npixRv 3 }
define npix1 $(npix1v[0])
define npix0 $(npix0v[0])
define npixR $(npixRv[0])
PNemax
#control code to find all of the max values of
pixels and save them in files?
#be sure rangeRD is set right
define rangeRD $($rang*$CD2_2)
PNe1max
PNe0max
PNeRmax
!rm "$!v.max.Im$!imag.R$!rang.max"
write
"$!v.max.Im$!imag.R$!rang.max" $npix1 $npix0 $npixR
PNe1max
#snippet code to count max # pixels possible for 1+
PNe regions
echo counting to max pixels for 1+ PNe
#calculate total number of pixels possible
# by looping through every pixel and seeing if
there are any PNe close enough to it
# that it could get counted, possibly, ever
define npix1 0
do px = 0, $($NAXIS1 - 1) {
do py = 0, $($NAXIS2 - 1) {
#check that it
is unmasked, if applicable
if ($(image($px,$py)) > -100) {
#convert to RA/DEC
define pDEC (($py - $CRPIX2 +
1)*$(ABS($CD2_2)) + $CRVAL2)
define pRA $($(-($px - $CRPIX1 + 1
))*$(ABS($CD1_1)/cosd($pDEC)) + $CRVAL1)
#check for any PNe (from RA/DEC vectors) within
2*rang+1 of pixel
set delRA = $pRA-RA
set delDEC = $pDEC-DEC
set counter =
((abs(delRA)<=$($($rangeRD*2)+$CD2_2)) &&
(abs(delDEC)<=$($($rangeRD*2)+$CD2_2)))? 1 : 0
if ($(sum(counter)) > 0) {
#found at least 1 PNe close
enough to this pixel
define npix1
$($npix1 + 1)
#end if some PNe found if statement
}
#end if pixel
is not masked
}
#end of y-pixels loop
}
echo $px / $($NAXIS1-1)
#end of x-pixels loop
}
echo max number of pixels with PNe is $npix1
PNe0max
#snippet of code to cound max # pixels possible for
0 PNe regions
echo counting to max pixels for 0 PNe
#loop through pixels, check if all PNe are far
enough away
define npix0 0
do px = 0, $($NAXIS1 - 1) {
do py = 0, $($NAXIS2 - 1) {
#check that it
is unmasked, if applicable
if ($(image($px,$py)) > -100) {
#convert to RA/DEC
define pDEC (($py - $CRPIX2 +
1)*$(ABS($CD2_2)) + $CRVAL2)
define pRA $($(-($px - $CRPIX1 + 1
))*$(ABS($CD1_1)/cosd($pDEC)) + $CRVAL1)
#check for any PNe (from RA/DEC
vectors) within rang of pixel
set delRA = $pRA-RA
set delDEC = $pDEC-DEC
set counter =
((abs(delRA)<=$rangeRD) &&
(abs(delDEC)<=$rangeRD))? 1 : 0
if ($(sum(counter)) == 0) {
#found NO PNe close enough to
this pixel
define npix0
$($npix0 + 1)
#end if some PNe found if statement
}
#end if pixel
is not masked
}
#end of y-pixels loop
}
echo $px / $($NAXIS1-1)
#end of x-pixels loop
}
echo max number of pixels without PNe is $npix0
PNeRmax
#ha, this is easy. max pixels in random boxes is
just image dimensions
#as long as pixels aren't masked
define npixR 0
do px = 0, $($NAXIS1 - 1) {
do py = 0, $($NAXIS2 - 1) {
if ($(image($px,$py)) > -100) {
define npixR $($npixR + 1)
}
}
#progress report
if ($($px/10) == int($px/10)) {
echo $px /
$($NAXIS1-1)
}
}
echo max number of pixels in random boxes is $npixR