Changeset 40693 for trunk

Timestamp:

Apr 13, 2019, 12:48:05 PM (7 years ago)

Author:

eugene

Message:

adding smear experiment

Location:

trunk/doc/release.2015/ps1.analysis/examples

Files:

• .mana (modified) (1 diff)
• moments.sh (modified) (2 diffs)

trunk/doc/release.2015/ps1.analysis/examples/.mana

@@ -84 +84 @@
 vstat e1
 vstat e2
+imsmooth.2d
+imsmooth.generic
+imsmooth
+fft2d
+input moments.sh
+input moments.sh
+opihi verbose on
+input moments.sh
+pwd
+!cat moments.sh
+input moments.sh
+opihi verbose off
+show.smear
+show.smear 1.0
+buffers
+input moments.sh
+show.smear 1.0
+input moments.sh
+show.smear 1.0
+buffers
+tv -n 0 -ch 1 f0 -0.01 0.5
+tv -n 0 -ch 1 f0 -0.01 1.0
+tv -n 0 -ch 2 psf -0.01 1.0
+tv -n 0 -ch 2 psf -0.01 0.5
+tv -n 0 -ch 2 fs -0.01 1.0
+stat fs
+tv -n 0 -ch 2 fs 0 2e-9
+tv im_R -0.01 0.1
+stat im_R
+tv im_R -0.0005 0.001
+tv im_I
+stat im_I
+tv kr_R -0.0005 0.001
+tv ff_R
+stat ff_R
+stat kr_R
+set ff_R = ff_R * 512^2
+set ff_I = ff_I * 512^2
+stat ff_R
+tv ff_R 0 0.05
+tv ff_I 0 0.05
+stat ff_I
+fft2d -inverse ff_R ff_I to ot_R ot_I
+stat ot_R
+set out_R = out_R * 512^2
+set out_R = ot_R * 512^2
+set out_I = ot_I * 512^2
+stat out_R
+stat out_I
+tv out_R
+set fs = fs * 512^4
+tv fs
+tv fs 0 152
+tv fs 0 512
+input moments.sh
+show.smear 1.0
+tv -n 0 -ch 2 fs 0 1
+tv -n 0 -ch 2 psf 0 1
+tv -n 0 -ch 2 psf 0 0.5
+tv -n 0 -ch 1 f0 0 0.5
+tv -n 0 -ch 3 fs 0 0.5
+star
+star f0 256 256 128
+star -q f0 256 256 128; $Mxx = ($SXg/2.355)^2; $Myy = ($SYg/2.355)^2;
+??
+star -q f0 256 256 128; $Mxx = ($SXg/2.355)^2; $Myy = ($SYg/2.355)^2; echo $Zcg $Mxx $SXYg $Myy
+star -q psf 256 256 128; $Mxx = ($SXg/2.355)^2; $Myy = ($SYg/2.355)^2; echo $Zcg $Mxx $SXYg $Myy
+star -q fs 256 256 128; $Mxx = ($SXg/2.355)^2; $Myy = ($SYg/2.355)^2; echo $Zcg $Mxx $SXYg $Myy
+echo {1256.63706807*25.1327415272}
+input moments.sh
+show.smear
+show.smear 1
+input moments.sh
+show.smear
+show.smear 1
+input moments.sh
+show.smear 1
+lim -n 2 theta e1_o; clear; box; plot theta e1_o -c blue; plot theta e1_s -c red
+vectors
+create t 0 720 5
+vectors
+create t 0 360 5
+vectors
+lim -n 2 t e1_o; clear; box; plot t e1_o -c blue; plot t e1_s -c red
+lim -n 2 t e1_o; clear; box; plot t e1_o -c blue; plot t e1_s -c red; plot t e2_o -c blue -pt box; plot t e2_s -c red -pt box
+lim -n 2 t e1_o; clear; box; plot t e1_o -c blue; plot t e1_s -c red; plot t e2_o -c blue -pt obox; plot t e2_s -c red -pt obox
+input moments.sh
+show.smear 1
+star psf {psf[][0]/2} {psf[0][]/2} 128
+$Mxx_p = ($SXg/2.355)^2; $Myy_p = ($SYg/2.355)^2; $Mxy_p = $SXYg
+lim -n 2 theta Mxx_s; clear; box; plot t Mxx_o -c blue; plot t Mxx_s -c red
+set Mxx_S = Mxx_o + $Mxx_p
+lim -n 2 theta Myy_s; clear; box; plot theta Myy_o -c blue; plot theta Myy_s -c red
+lim -n 2 theta 0 200; clear; box; plot theta Myy_o -c blue; plot theta Myy_s -c red
+lim -n 2 theta 0 150; clear; box; plot theta Myy_o -c blue; plot theta Myy_s -c red
+set dMxx = Mxx_s - Mxx_o
+set dMxy = Mxy_s - Mxy_o
+set dMyy = Myy_s - Myy_o
+lim -n 2 theta dMxx; clear; box; plot theta dMxx
+line -c red $Mxx_p 0 to $Mxx_p 360
+line -c red 0 $Mxx_p to 360 $Mxx_p
+lim -n 2 theta dMxx; clear; box; plot theta dMxx; line -c red 0 $Mxx_p to 360 $Mxx_p
+lim -n 2 theta dMyy; clear; box; plot theta dMyy; line -c red 0 $Myy_p to 360 $Myy_p
+vstat dMyy
+lim -n 2 theta 35 37; clear; box; plot theta dMyy; line -c red 0 $Myy_p to 360 $Myy_p
+lim -n 2 theta dMxx; clear; box; plot theta dMxx; line -c red 0 $Mxx_p to 360 $Mxx_p
+lim -n 2 theta dMxy; clear; box; plot theta dMxy; line -c red 0 $Mxy_p to 360 $Mxy_p
+lim -n 1 theta Mxy_s; clear; box; plot theta Mxy_s -c red -pt cir -sz 2; plot theta Mxy_o -c blue -pt cir -sz 2; line -c red 0 $Mxy_p to 360 $Mxy_p
+input moments.sh
+show.smear 1
+lim -n 1 theta Mxy_s; clear; box; plot theta Mxy_s -c red -pt cir -sz 2; plot theta Mxy_o -c blue -pt cir -sz 2; line -c red 0 $Mxy_p to 360 $Mxy_p
+set dMxy = Mxy_s - Mxy_o
+lim -n 2 theta dMxy; clear; box; plot theta dMxy; line -c red 0 $Mxy_p to 360 $Mxy_p

trunk/doc/release.2015/ps1.analysis/examples/moments.sh

@@ +1 @@
+
+macro imconvolve
+  if ($0 != 4) 
+    echo "USAGE: imconvolve (image) (kernel) (output)"
+    break
+  end
+
+  # if kernel is not the same size as input, pad
+  if ($1[][0] != $2[][0]) 
+    echo "image sizes do not match"
+    break
+  end
+  if ($1[0][] != $2[0][]) 
+    echo "image sizes do not match"
+    break
+  end
+
+  local Nx Ny
+  $Nx = $1[][0]
+  $Ny = $1[0][]
+
+  fft2d $1 0 to im_R im_I
+  fft2d $2 0 to kr_R kr_I
+
+  $Npix = $Nx*$Ny
+  set ff_R = $Npix*(im_R * kr_R - im_I * kr_I)
+  set ff_I = $Npix*(im_I * kr_R + im_R * kr_I)
+
+  fft2d -inverse ff_R ff_I to ot_R ot_I
+  set tmppow = $Npix*sqrt(ot_R^2 + ot_I^2)
+
+  extract tmppow $3       0       0 {$Nx/2} {$Ny/2} {$Nx/2} {$Ny/2} $Nx $Ny
+  extract tmppow $3 {$Nx/2}       0 {$Nx/2} {$Ny/2}       0 {$Ny/2} $Nx $Ny
+  extract tmppow $3       0 {$Ny/2} {$Nx/2} {$Ny/2} {$Nx/2}       0 $Nx $Ny
+  extract tmppow $3 {$Nx/2} {$Ny/2} {$Nx/2} {$Ny/2}       0       0 $Nx $Ny
+end
 
 macro show.pol
@@ -30 +66 @@
 end
 
+macro show.smear
+  if ($0 != 2) 
+    echo "USAGE: show.smear (sigma)"
+    break
+  end
+
+  mcreate z 512 512
+  set x = xramp(z) - z[][0]/2
+  set y = yramp(z) - z[0][]/2
+
+  set f0 = 10*exp(-0.5*((x/2)^2 + (y/10)^2))
+  set psf_r = exp(-0.5*((x/2)^2 + (y/6)^2))
+  rotate psf_r 30.0
+  extract psf_r psf {psf_r[][0]/2 - 256} {psf_r[0][]/2 - 256} 512 512 0 0 512 512
+  star -q psf {psf[][0]/2} {psf[0][]/2} 128
+  $Mxx_p = ($SXg/2.355)^2
+  $Myy_p = ($SYg/2.355)^2
+  $Mxy_p = $SXYg
+
+  delete -q e1_s e2_s e1_o e2_o theta Mxx_s Mxy_s Myy_s Mxx_o Mxy_o Myy_o
+  for rot 0 360 5
+    set fr = f0
+    rotate fr $rot
+    delete -q frs
+    extract fr frs {fr[][0]/2 - 256} {fr[0][]/2 - 256} 512 512 0 0 512 512
+
+    imconvolve frs psf fs
+
+    star -q fs {fs[][0]/2} {fs[0][]/2} 128
+
+    $Mxx = ($SXg/2.355)^2
+    $Myy = ($SYg/2.355)^2
+    $Mxy = $SXYg
+
+    $T = $Mxx + $Myy
+    $e1 = ($Mxx - $Myy) / $T
+    $e2 = 2*$Mxy / $T
+
+    concat $rot theta
+    concat $e1 e1_s
+    concat $e2 e2_s
+    concat $Mxx Mxx_s
+    concat $Mxy Mxy_s
+    concat $Myy Myy_s
+
+    # fr is the rotate version of f0
+    star -q fr {fr[][0]/2} {fr[0][]/2} 128
+
+    $Mxx = ($SXg/2.355)^2
+    $Myy = ($SYg/2.355)^2
+    $Mxy = $SXYg
+
+    $T = $Mxx + $Myy
+    $e1 = ($Mxx - $Myy) / $T
+    $e2 = 2*$Mxy / $T
+
+    concat $e1 e1_o
+    concat $e2 e2_o
+    concat $Mxx Mxx_o
+    concat $Mxy Mxy_o
+    concat $Myy Myy_o
+  end
+
+  lim -n graph -1.5 1.5 -1.5 1.5; clear; box; plot -pt cir -sz 2 -c blue e1_o e2_o; plot -pt cir -sz 2 -c red e1_s e2_s; 
+end
 
 macro test.rot