Skip to content

Instantly share code, notes, and snippets.

@jberger
Created June 28, 2012 17:43
Show Gist options
  • Save jberger/3012794 to your computer and use it in GitHub Desktop.
Save jberger/3012794 to your computer and use it in GitHub Desktop.
PGPS gaussian width app
#!/usr/bin/env perl
use strict;
use warnings;
use PDL;
use PDL::Image2D;
use PDL::Fit::Gaussian;
use Prima;
use PDL::Graphics::Prima::Simple -sequential;
use File::chdir;
use List::Util;
use constant {
pi => 4*atan2(1,1),
wconv => 1 / ( 2 * sqrt( log(2) ) ),
};
my $folder = shift;
local $CWD = $folder;
my $cache_name = 'data';
my $tag = 'fourier';
my $nd_tag = qr/g(\d+)/;
my $handle_backgrounds = 0;
my $too_large = 900;
my $laser = sub{ sin(pi / 180 * 2 * shift)**2 };
my $files;
if ( not -e $cache_name ) {
$files = analyze_folder();
# cache data
require Data::Dumper;
open my $fh, '>', $cache_name;
print $fh Data::Dumper::Dumper($files);
} else {
$files = do $cache_name;
}
# print for humans
use Data::Printer;
p $files;
line_plot(make_piddles($files));
sub analyze_folder {
# find all fit(s) files
my %files = map { $_, { file => $_ } } do {
opendir( my $dh, $CWD );
grep { /\.fits?/ }
readdir $dh;
};
# keep only certain tags (e.g. fourier)
if ($tag) {
my @wrong_tag = grep { not /$tag/i } keys %files;
delete @files{@wrong_tag};
}
# remove background images
unless ($handle_backgrounds) {
my @bg = grep { /bg/ } keys %files;
delete @files{@bg};
}
# parse filenames
foreach (keys %files) {
# polarizer info
if ( /(m?)(\d{2,})/ ) {
$files{$_}{polarizer} = $2 * ($1 ? -1 : 1); # m indicates negative number
$files{$_}{laser} = $laser->($files{$_}{polarizer});
}
# nd info
if ( $_ =~ $nd_tag ) {
$files{$_}{nd} = $1;
} else {
$files{$_}{nd} = 0;
}
}
foreach (sort {$files{$a}{laser} <=> $files{$b}{laser}} keys %files) {
my @w = analyze($_);
$files{$_}{w} = \@w;
$files{$_}{w_av} = 0.5 * List::Util::sum @w;;
}
return \%files;
}
sub make_piddles {
my $files = shift;
my (@laser, @w_av);
foreach (sort {$files->{$a}{laser} <=> $files->{$b}{laser}} keys %$files) {
next if $files->{$_}{w_av} > $too_large;
push @laser, $files->{$_}{laser};
push @w_av, $files->{$_}{w_av};
}
my $pdl_laser = pdl \@laser;
my $pdl_w_av = pdl \@w_av;
return $pdl_laser, $pdl_w_av;
}
sub analyze {
my $file = shift;
return unless $file =~ /\.fits?$/;
print "processing file: $file\n";
my $pdl = rfits $file;
my ($x_dim, $y_dim) = $pdl->dims;
my ($val, $max_x, $max_y) = $pdl->max2d_ind;
my ($wx, $wy);
plot(
-data => ds::Grid(
$pdl,
x_bounds => [1,$x_dim],
y_bounds => [1,$y_dim],
plotType => pgrid::Matrix(
palette => pal::BlackToWhite,
),
),
-max => ds::Pair(
$max_x, $max_y,
plotTypes => ppair::Crosses,
colors => pdl(255, 0, 0)->rgb_to_color,
),
onMouseClick => sub {
my ($self, $button, undef, $x, $y) = @_;
my $left = $button & mb::Left;
my $right = $button & mb::Right;
return unless ($left or $right);
if ( $left ) {
$x = int $self->x->pixels_to_reals($x);
$y = int $self->y->pixels_to_reals($y); #-- #highlight fix
} else {
($x, $y) = ($max_x, $max_y);
}
($wx, $wy) = xy_fit($pdl, $x, $y);
$self->get_parent->close;
},
);
return ($wx, $wy);
}
sub xy_fit {
my ($pdl, $x, $y) = @_;
my $x_lineout = $pdl->slice(",($y)");
my $y_lineout = $pdl->slice("($x),");
my (undef, undef, $fwhm_x) = fitgauss1d(
$x_lineout->sequence * 5.4,
$x_lineout,
);
my (undef, undef, $fwhm_y) = fitgauss1d(
$y_lineout->sequence * 5.4,
$y_lineout,
);
return
map { wconv * $_ }
map {eval{$_->isa('PDL')} ? $_->list : $_}
($fwhm_x, $fwhm_y);
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment