/packmule_raid/steven/PNe12/PNeError.sm
PNeError
#two part program
# -first to compute the values for the
values for the error bars
# -second to graph
# -the average values
for each image/range over all iterations
# -the error bars!!
#first set up vectors of all possibilities
set im = { 3.fits 3m.fits 4.fits 4m.fits 7.fits
7m.fits Sub.fits Subm.fits }
data v.v
read { vnum 1 }
define v $(vnum[0])
data i.i
read { inum 1 }
define itt $(inum[0])
set it = 1,$itt
#set ra = 10,60,5
set ra = { 1 3 5 7 9 11 13 15 17 20 25 30 35 40 45
50 55 60 }
#temp for testing
# set im = {3.fits}
# set ra = 10,60,5
#separate main runings, first by imagee
foreach imag im {
#prep vectors for graphing
set dimen(wGm) = 0
set dimen(wGs) = 0
set dimen(wGu) = 0
set dimen(wGl) = 0
set dimen(woGm) = 0
set dimen(woGs) = 0
set dimen(woGu) = 0
set dimen(woGl) = 0
set dimen(randGm) = 0
set dimen(randGs) = 0
set dimen(randGu) = 0
set dimen(randGl) = 0
set dimen(wMt) = 0
set dimen(woMt) = 0
set dimen(randMt) = 0
#add on for delta support
set dimen(DwGm) = 0
set dimen(DwGs) = 0
set dimen(DwGu) = 0
set dimen(DwGl) = 0
set dimen(DwoGm) = 0
set dimen(DwoGs) = 0
set dimen(DwoGu) = 0
set dimen(DwoGl) = 0
set dimen(DwwoGm) = 0
set dimen(DwwoGs) = 0
set dimen(DwwoGu) = 0
set dimen(DwwoGl) = 0
#next loop is by image
foreach rang ra {
#call
computational program - it loops through ranges
PNeErrorc
#end of ranges loop
}
#call graphing program - now all
ranges for imag have been done
PNeErrorg
#don't ask why but it's happier
if you run it twice
#as time permits, leave the next
line un-commented-out
PNeErrorg
#end of images loop
}
echo converting to gifs/uploading to the web and
Virgo/Web folder...
qc
echo done
PNeErrorc
#this program computes error bars for the PNe7 run,
based on .med files
# which should be labeled 7.#.Im3.fits.R10.med
# where # is iteration, 1 to ____
# 3.fits file name,
3,3m,4,4m,7,7m,Sub,Subm
# R10 is range, 10 to 60, by 5s
# file will look like:
# write
"$!v.$!itn.Im$!imag.R$!rang.med" "$!imag
$!rang $!wM $!woM
$!randM "
set DIMEN(wMt)
= 0
set
DIMEN(woMt) = 0
set
DIMEN(randMt) = 0
#innermost
loop is by iteration number
foreach itn it
{
#clear variables
set dimen(wMv) = 0
set dimen(woMv) = 0
set dimen(randMv) = 0
#read in data IF file exists
if (is_file($v.$itn.Im$imag.R$rang.med)) {
data $v.$itn.Im$imag.R$rang.med
read { imchk 1.s rachk 2 wMv 3
woMv 4 randMv 5 }
#store values in temp vectors
set wMt = wMt concat $(wMv[0])
set woMt = woMt concat $(woMv[0])
set randMt = randMt concat
$(randMv[0])
}
#end of
iterations loop
}
#all
iterations for one range/one image - stats time
#compute
deltas for this range
set Dw = wMt -
randMt
set Dwo = woMt
- randMt
set Dwwo = wMt
- woMt
#finding sigma:
stats wMt
wMean wSigma wK
stats woMt
woMean woSigma woK
stats randMt
randMean randSigma randK
stats Dw
DwMean DwSigma DwK
stats Dwo
DwoMean DwoSigma DwoK
stats Dwwo
DwwoMean DwwoSigma DwwoK
#stats_med
does not output quartiles so copy/do manually
sort {wMt}
define wMed
(wMt[DIMEN(wMt)/2])
define wLQ
(wMt[0.25*DIMEN(wMt)])
define wUQ
(wMt[0.75*DIMEN(wMt)])
sort {woMt}
define woMed
(woMt[DIMEN(woMt)/2])
define woLQ
(woMt[0.25*DIMEN(woMt)])
define woUQ
(woMt[0.75*DIMEN(woMt)])
sort {randMt}
define randMed
(randMt[DIMEN(randMt)/2])
define randLQ
(randMt[0.25*DIMEN(randMt)])
define randUQ
(randMt[0.75*DIMEN(randMt)])
sort {Dw}
define DwMed
(Dw[DIMEN(Dw)/2])
define DwLQ
(Dw[0.25*DIMEN(Dw)])
define DwUQ
(Dw[0.75*DIMEN(Dw)])
sort {Dwo}
define DwoMed
(Dwo[DIMEN(Dwo)/2])
define DwoLQ
(Dwo[0.25*DIMEN(Dwo)])
define DwoUQ
(Dwo[0.75*DIMEN(Dwo)])
sort {Dwwo}
define DwwoMed
(Dwwo[DIMEN(Dwwo)/2])
define DwwoLQ
(Dwwo[0.25*DIMEN(Dwwo)])
define DwwoUQ
(Dwwo[0.75*DIMEN(Dwwo)])
#write results
to data file
!rm
"$!v.0.Im$!imag.R$!rang.stats"
#writerand1Gm
header
write
"$!v.0.Im$!imag.R$!rang.stats" "#wPNe woPNe rand Dw Dwo Dwoo - data
rows go: (1)Mean, (2)Sigma, (3)Kurtosis, (4)Median, (5)Lower Quartile,
(6)Upper Quartile"
write
"$!v.0.Im$!imag.R$!rang.stats" $wMean $woMean $randMean $DwMean
$DwoMean $DwwoMean
write
"$!v.0.Im$!imag.R$!rang.stats" $wSigma $woSigma $randSigma $DwSigma
$DwoSigma $DwwoSigma
write
"$!v.0.Im$!imag.R$!rang.stats" $wK $woK $randK $DwK $DwoK $DwwoK
write
"$!v.0.Im$!imag.R$!rang.stats" $wMed $woMed $randMed $DwMed $DwoMed
$DwwoMed
write
"$!v.0.Im$!imag.R$!rang.stats" $wLQ $woLQ $randLQ $DwLQ $DwoLQ $DwwoLQ
write
"$!v.0.Im$!imag.R$!rang.stats" $wUQ $woUQ $randUQ $DwUQ $DwoUQ $DwwoUQ
#store
variables for graphing soon
set wGm = wGm
CONCAT $wMean
set wGs = wGs
CONCAT $wSigma
set wGu = wGu
CONCAT $wUQ
set wGl = wGl
CONCAT $wLQ
set woGm =
woGm CONCAT $woMean
set woGs =
woGs CONCAT $woSigma
set woGu =
woGu CONCAT $woUQ
set woGl =
woGl CONCAT $woLQ
set randGm =
randGm CONCAT $randMean
set randGs =
randGs CONCAT $randSigma
set randGu =
randGu CONCAT $randUQ
set randGl =
randGl CONCAT $randLQ
#add'l delta
support
set DwGm =
DwGm CONCAT $DwMean
set DwGs =
DwGs CONCAT $DwSigma
set DwGu =
DwGu CONCAT $DwUQ
set DwGl =
DwGl CONCAT $DwLQ
set DwoGm =
DwoGm CONCAT $DwoMean
set DwoGs =
DwoGs CONCAT $DwoSigma
set DwoGu =
DwoGu CONCAT $DwoUQ
set DwoGl =
DwoGl CONCAT $DwoLQ
set DwwoGm =
DwwoGm CONCAT $DwwoMean
set DwwoGs =
DwwoGs CONCAT $DwwoSigma
set DwwoGu =
DwwoGu CONCAT $DwwoUQ
set DwwoGl =
DwwoGl CONCAT $DwwoLQ
PNeErrorg
expand 1.0001
define TeX_strings 1
#graphin' time (tm)
#this is run once for each image (3,3m,4,4m,etc)
#can now use vectors defined above in PNeErrorc to
graph things
#data in vectors is one row / range value, ranges
are 10-60, by 5s
#will graph in terms of arcseconds, not pixels
set rangesG = ra*1.45
device postencap '$!v.0.med$imag.we.ps'
# device x11
erase
#top window left: 1+ PNe Medians
window 2 -3 1 3
ctype black
vecminmax wGm wmin wmax
lim -1 90 $($wmin-$($wmin*.5))
$($wmax+$($wmax*.5))
box
connect rangesG wGm
# #temp for initial testing
# relocate $(rangesG[0]) $(wGm[0])
# #dot $($wmin-$($wmin*.5))
$($wmax+$($wmax*.5))
ylabel Random w/ 1+ PNe
#sigma error bars
ctype blue
set sig = abs(wGs)
errorbar rangesG wGm sig 2
errorbar rangesG wGm sig 4
#U/L quartile error bars
ctype red
set R2u = abs(wGu - wGm)
set R2l = abs(wGl - wGm)
errorbar rangesG wGm R2u 2
errorbar rangesG wGm R2l 4
#basic labels
lim 10 90 0 10.2
ctype black
relocate 25 11.6
label $imag -
ctype red
relocate 110 11.6
label red-upper/lower quartile
ctype blue
relocate 70 11.6
label blue-sigma
ctype black
relocate 25 10.5
label average medians
relocate 143 10.5
label delta plots
relocate 87 10.3
label \it{$(dimen(randMt))+ trials}
#top middle window: Random PNe Medians
window 2 -3 1 2
ctype black
vecminmax randGm randmin randmax
lim -1 90
$($randmin-$($randmin*.5)) $($randmax+$($randmax*.5))
box
connect rangesG randGm
# #temp for initial testing
# relocate $(rangesG[0])
$(randGm[0])
# dot
ylabel Random
#sigma error bars
ctype blue
set sig = abs(randGs)
errorbar rangesG randGm sig 2
errorbar rangesG randGm sig 4
#U/L quartile error bars
ctype red
set Ru = abs(randGu - randGm)
set Rl = abs(randGl - randGm)
errorbar rangesG randGm Ru 2
errorbar rangesG randGm Rl 4
#bottom middle window: 0 PNe Medians
window 2 -3 1 1
ctype black
vecminmax woGm womin womax
lim -1 90 $($womin-$($womin*.5))
$($womax+$($womax*.5))
box
connect rangesG woGm
# #temp for initial testing
# relocate $(rangesG[0]) $(woGm[0])
# dot
ylabel Random w/ 0 PNe
xlabel half-box width (arcseconds)
#sigma error bars
ctype blue
set sig = abs(woGs)
errorbar rangesG woGm sig 2
errorbar rangesG woGm sig 4
#U/L quartile error bars
ctype red
set wou = abs(woGu-woGm)
set wol = abs(woGl-woGm)
errorbar rangesG woGm wou 2
errorbar rangesG woGm wol 4
#do the deltas!
##calculate deltas - ANTIQUATED
##set Rw = wGm - randGm
##set Rwo = woGm - randGm
##set Rwwo = wGm - woGm
#find limits
vecminmax DwGm wmin wmax
vecminmax DwoGm womin womax
vecminmax DwwoGm wwomin wwomax
set dimen(chk) = 6
set chk[0]=$wmin
set chk[1]=$wmax
set chk[2]=$womin
set chk[3]=$womax
set chk[4]=$wwomin
set chk[5]=$wwomax
vecminmax chk ymin ymax
#if ($twomin < $onmin && $twomin <
$offmin) { define ymin $twomin }
#if ($onmin < $twomin && $onmin <
$offmin) { define ymin $onmin }
#if ($offmin < $twomin && $offmin <
$onmin) { define ymin $offmin }
#if ($twomax > $onmax && $twomax >
$offmax) { define ymax $twomax }
#if ($onmax > $twomax && $onmax >
$offmax) { define ymax $onmax }
#if ($offmax > $twomax && $offmax >
$onmax) { define ymax $offmax }
#add some buffer around the edges for error bars
define ymax $($ymax + $(abs($ymax-$ymin)/2))
define ymin $($ymin - $(abs($ymax-$ymin)/2))
#deltas on the right side
#top window : delta (R - 1)
window 2 -3 2 3
ctype black
lim -1 90 $ymin $ymax
box
connect rangesG DwGm
ylabel Delta(1 - R)
#sigma error bars
ctype blue
set sig = abs(DwGs)
errorbar rangesG DwGm sig 2
errorbar rangesG DwGm sig 4
#U/L quartile error bars
ctype red
set u = abs(DwGu-DwGm)
set l = abs(DwGl-DwGm)
errorbar rangesG DwGm u 2
errorbar rangesG DwGm l 4
#middle window : delta (R - 0)
window 2 -3 2 2
ctype black
lim -1 90 $ymin $ymax
box
connect rangesG DwoGm
ylabel Delta(0 - R)
#sigma error bars
ctype blue
set sig = abs(DwoGs)
errorbar rangesG DwoGm sig 2
errorbar rangesG DwoGm sig 4
#U/L quartile error bars
ctype red
set u = abs(DwoGu - DwoGm)
set l = abs(DwoGl - DwoGm)
errorbar rangesG DwoGm u 2
errorbar rangesG DwoGm l 4
#bottom window : delta (1 - 0)
window 2 -3 2 1
ctype black
lim -1 90 $ymin $ymax
box
connect rangesG DwwoGm
ylabel Delta(1 - 0)
xlabel half-box width (arcseconds)
#sigma error bars
ctype blue
set sig = abs(DwwoGs)
errorbar rangesG DwwoGm sig 2
errorbar rangesG DwwoGm sig 4
#U/L quartile error bars
ctype red
set u = abs(DwwoGu - DwwoGm)
set l = abs(DwwoGl - DwwoGm)
errorbar rangesG DwwoGm u 2
errorbar rangesG DwwoGm l 4
qc
!convert 12.0.med3.fits.we.ps
12.0.med3.fits.we.ps.gif
!convert 12.0.med3m.fits.we.ps
12.0.med3m.fits.we.ps.gif
!convert 12.0.med4.fits.we.ps
12.0.med4.fits.we.ps.gif
!convert 12.0.med4m.fits.we.ps
12.0.med4m.fits.we.ps.gif
!convert 12.0.med7.fits.we.ps
12.0.med7.fits.we.ps.gif
!convert 12.0.med7m.fits.we.ps
12.0.med7m.fits.we.ps.gif
!convert 12.0.medSub.fits.we.ps
12.0.medSub.fits.we.ps.gif
!convert 12.0.medSubm.fits.we.ps
12.0.medSubm.fits.we.ps.gif
!cp *.gif /astroweb_data/web/steven/PNe12/
!cp *.gif ~/Virgo/Web/PNe12/