Have you ever found yourself with an obnoxious blot in an otherwise beautiful sky? You want to remove the blot, but filling in the sky texture is not easy. The sky is not a simple color, and changes with the angle to the horizon.
Here is an example (this is a thumbnail; click on it to see the full picture):
The plane is ugly—the picture would look much better without it. It is easy enough to erase the plane:
However, if we just fill the erased area with a solid blue color, it will look patched. Here is a better solution. I have written a program, fill_erasures, which replaces the erased area with the texture from its edges. It does a statistical analysis of the picture, then fills the erased area by interpolation. Here is the result:
Here is how you use the program. Take your original picture, erase the unwanted sky object and convert it to PNG format. Run fill_erasures with the input file as the first argument and an output file as the second argument. Load the input file into an image editing program such as The Gimp, and load the output file as a layer underneath the input file. The interpolated sky texture will replace the erased parts of the input file.
fill_erasures is written in perl, so it should run on GNU/Linux, Apple OS/X and Microsoft Windows, since these all have perl available. However, I have only tested it on GNU/Linux
Here is the code. You should be able to cut and paste it into a file and execute it using your perl interpreter. If you prefer, you can fetch it from this URL: https://www.systemeyescomputerstore.com/texture_processing/fill_erasures . Try it on some of your own pictures and let me know how well it works for you. There is help text at the end.
#!/usr/local/bin/perl
use 5.010;
use GD;
use Math::Random;
use Getopt::Long;
use Pod::Usage;
my $man = 0;
my $help = 0;
my $patch_size = 10;
my $blend_patches = 0;
#
# Parse options.
#
my $getopt_result = GetOptions (
"help|?" => \$help, "man" => \$man,
"patch_size=i" => \$patch_size,
"blend_patches" => \$blend_patches,
"verbose+" => \$verbose) or pod2usage(2);
pod2usage(1) if $help;
pod2usage(-exitstatus => 0, -verbose => 2) if $man;
#
# Accept parameters.
#
my $in_name = shift (@ARGV);
my $out_name = shift (@ARGV);
if (!defined(in_name) or $in_name eq "") {
die "input file name is required; use --help for help.\n";
}
if (!defined(out_name) or $out_name eq "") {
die "output file name is required; use --help for help.\n";
}
if (! (-r $in_name)) {
die "input file ", $in_name, " cannot be read.\n";
}
my $in_image = GD::Image->newFromPng($in_name, 1);
my $image_width = $in_image->width;
my $image_height = $in_image->height;
if ($verbose > 0) {
print "Image ", $in_name, " is ", $image_width, " by ", $image_height, ".\n";
print "Patch size is ", $patch_size, " by ", $patch_size, " pixels.\n";
}
#
# The GD package doesn't provide a way to extract the alpha
# value of a pixel, so we cheat by knowing the internal
# representation of an index. Any pixels off the image
# are transparent.
#
sub getRGBA {
my ($x, $y) = @_;
if (($x < 0) or ($x >= $image_width) or ($y < 0) or ($y >= $image_height)) {
if ($verbose > 3) {
print "[", $x, ", ", $y, "] = off screen.\n";
}
return (0, 0, 0, 127);
}
my $index = $in_image->getPixel ($x, $y);
my ($red_val, $green_val, $blue_val) = $in_image->rgb($index);
my $alpha_val = int($index / (2**24));
if ($verbose > 3) {
print "[", $x, ", ", $y, "] = (", $red_val, ", ", $green_val, ", ",
$blue_val, ", ", $alpha_val, ");\n";
}
return ($red_val, $green_val, $blue_val, $alpha_val);
}
#
# Compute the statistics for each patch of the picture.
# Patches with less than patch_size opaque pixels are considered
# erased. Thus we don't lose the last row or column if the picture
# width or height is one more than a multiple of patch_size.
#
# Internally we use scRGB to represent colors.
# This is a linear color system based on sRGB, from IEC 61966-2-2.
#
my %original_patches;
my %replacement_patches;
for (my $outer_x = 0; $outer_x < $image_width; $outer_x = $outer_x + $patch_size) {
for (my $outer_y = 0; $outer_y < $image_height; $outer_y = $outer_y + $patch_size) {
my $pixel_count = 0;
my $red_sum = 0;
my $green_sum = 0;
my $blue_sum = 0;
for (my $x = $outer_x; $x < $outer_x + $patch_size; ++$x) {
for (my $y = $outer_y; $y < $outer_y + $patch_size; ++$y) {
my ($red_val, $green_val, $blue_val, $alpha_val) = &getRGBA($x, $y);
if ($alpha_val == 0) {
$pixel_count = $pixel_count + 1;;
$red_val = &to_scrgb ($red_val);
$green_val = &to_scrgb ($green_val);
$blue_val = &to_scrgb ($blue_val);
$red_sum += $red_val;
$green_sum += $green_val;
$blue_sum += $blue_val;
}
}
}
if ($pixel_count >= $patch_size) {
my $red_mean = $red_sum / $pixel_count;
my $green_mean = $green_sum / $pixel_count;
my $blue_mean = $blue_sum / $pixel_count;
my $red_dev = 0;
my $green_dev = 0;
my $blue_dev = 0;
for (my $y = $outer_y; $y < $outer_y + $patch_size; $y++) {
for (my $x = $outer_x; $x < $outer_x + $patch_size; $x++) {
my ($red_val, $green_val, $blue_val, $alpha_val) = &getRGBA ($x, $y);
if ($alpha_val eq 0) {
$red_val = &to_scrgb ($red_val);
$green_val = &to_scrgb ($green_val);
$blue_val = &to_scrgb ($blue_val);
my $red_deviation = ($red_mean - $red_val);
$red_dev += ($red_deviation)**2;
my $green_deviation = ($green_mean - $green_val);
$green_dev += ($green_deviation)**2;
my $blue_deviation = ($blue_mean - $blue_val);
$blue_dev += ($blue_deviation)**2;
}
}
}
my $red_stddev = sqrt($red_dev / $pixel_count);
my $green_stddev = sqrt($green_dev / $pixel_count);
my $blue_stddev = sqrt($blue_dev / $pixel_count);
$original_patches{"${outer_x} ${outer_y}"} =
"${red_mean} ${green_mean} ${blue_mean} " .
"${red_stddev} ${green_stddev} ${blue_stddev}";
if ($verbose > 1) {
print "[", $outer_x, ", ", $outer_y, "] is (", $red_mean, ", ",
$green_mean, ", ", $blue_mean, "), ",
$red_stddev, ", ", $green_stddev, ", ", $blue_stddev, "\n";
}
} else {
$replacement_patches{"${outer_x} ${outer_y}"} = "0 0 0 0 0 0 0";
}
}
}
if ($verbose > 0) {
my $replacement_patch_count = keys %replacement_patches;
print "There are ", $replacement_patch_count, " erased patches to be filled in.\n";
}
#
# Go through the patches, filling in the erased patches by averaging
# with neighbors.
#
sub find_neighbors {
#
# Return a hash of the up to four possible neighbors of this patch.
# We find the nearest original patch in four directions. This yields
# up to four neighbors.
#
my ($x, $y) = @_;
my %neighbors;
if ($verbose > 3) {
print "Neighbors of [", $x, ", ", $y, "].\n";
}
foreach (keys %original_patches) {
my ($neighbor_x, $neighbor_y) = split " ", $_;
if ($verbose > 3) {
print " considering [", $neighbor_x, ", ", $neighbor_y, "]: ";
}
my $delta_x = $neighbor_x - $x;
my $delta_y = $neighbor_y - $y;
my $distance = sqrt (($delta_x**2) + ($delta_y**2));
if ($verbose > 3) {
print $distance, " ", $delta_x, " ", $delta_y, "\n";
}
if (($delta_x >= 0) and ($delta_y < 0)) {
if ($verbose > 3) {
print " above+right\n";
}
if (exists($neighbors{"above+right"})) {
my ($oldn_distance, $oldn_x, $oldn_y) = split " ", $neighbors{"above+right"};
if ($distance < $oldn_distance) {
$neighbors{"above+right"} = $distance . " " . $neighbor_x . " " . $neighbor_y;
}
} else {
$neighbors{"above+right"} = $distance . " " . $neighbor_x . " " . $neighbor_y;
}
}
if (($delta_x >= 0) and ($delta_y >= 0)) {
if ($verbose > 3) {
print " below+right\n";
}
if (exists($neighbors{"below+right"})) {
my ($oldn_distance, $oldn_x, $oldn_y) = split " ", $neighbors{"below+right"};
if ($distance < $oldn_distance) {
$neighbors{"below+right"} = $distance . " " . $neighbor_x . " " . $neighbor_y;
}
} else {
$neighbors{"below+right"} = $distance . " " . $neighbor_x . " " . $neighbor_y;
}
}
if (($delta_x < 0) and ($delta_y >= 0)) {
if ($verbose > 3) {
print " below+left\n";
}
if (exists($neighbors{"below+left"})) {
my ($oldn_distance, $oldn_x, $oldn_y) = split " ", $neighbors{"below+left"};
if ($distance < $oldn_distance) {
$neighbors{"below+left"} = $distance . " " . $neighbor_x . " " . $neighbor_y;
}
} else {
$neighbors{"below+left"} = $distance . " " . $neighbor_x . " " . $neighbor_y;
}
}
if (($delta_x < 0) and ($delta_y < 0)) {
if ($verbose > 3) {
print " above+left\n";
}
if (exists($neighbors{"above+left"})) {
my ($oldn_distance, $oldn_x, $oldn_y) = split " ", $neighbors{"above+left"};
if ($distance < $oldn_distance) {
$neighbors{"above+left"} = $distance . " " . $neighbor_x . " " . $neighbor_y;
}
} else {
$neighbors{"above+left"} = $distance . " " . $neighbor_x . " " . $neighbor_y;
}
}
}
if ($verbose > 3) {
while ( my ($key, $val) = each %neighbors) {
print " ", $key, " => ", $val, "\n";
}
}
return %neighbors;
}
#
# For each patch that needs replacement, find its surrounding neighbors and
# fill in its statistics based on their distance.
#
my %all_patches = %original_patches;
while ( my ($position, $statistics) = each %replacement_patches) {
my ($x, $y) = split " ", $position;
my %neighbors = &find_neighbors ($x, $y);
my @statistics = (0,0,0,0,0,0);
my $weight_sum = 0;
while ( my ($direction, $distance_and_position) = each %neighbors) {
my ($distance, $neighbor_x, $neighbor_y) = split " ", $distance_and_position;
my $weight = (1.0 / $distance);
my @neighbor_statistics = split " ", $original_patches{"${neighbor_x} ${neighbor_y}"};
if ($verbose > 2) {
print " [", $neighbor_x, ", ", $neighbor_y, "] --> ",
join " ", @neighbor_statistics, " ", $direction, " at ", $distance, "\n";
}
for (my $stat = 0; $stat < 6; $stat++) {
@statistics[$stat] = @statistics[$stat] + (@neighbor_statistics[$stat] * $weight);
}
$weight_sum = $weight_sum + $weight;
}
for (my $stat = 0; $stat < 6; $stat++) {
@statistics[$stat] = @statistics[$stat] / $weight_sum;
}
$all_patches{"${x} ${y}"} = join " ", @statistics;
if ($verbose > 1) {
print "[", $x, ", ", $y, "] gets ", $all_patches{"${x} ${y}"}, "\n";
}
}
#
# If requested, for each pixel in the output image, determine its statistics
# by averaging the statistics of the nearby patches, weighted
# by their distance. Otherwise, a pixel just gets its patch's statistics.
#
my $out_image = GD::Image->new($image_width, $image_height,1);
if ($blend_patches == 1) {
if ($verbose > 0) {
print "Blending patches.\n";
}
my $max_distance = sqrt(2) * sqrt ((($patch_size / 2)**2) + (($patch_size / 2) ** 2));
my %pixel_stats;
if ($verbose > 2) {
print "Max distance = ", $max_distance, "\n";
}
for ($x = 0; $x < $image_width; $x++) {
PIXEL:
for ($y = 0; $y < $image_height; $y++) {
if ($verbose > 2) {
print "Starting [", $x, ", ", $y, "]\n";
}
my %nearby_patches;
foreach (keys %all_patches) {
my ($patch_x, $patch_y) = split " ", $_;
$patch_x = $patch_x + ($patch_size / 2);
$patch_y = $patch_y + ($patch_size / 2);
my $distance = sqrt ((($x - $patch_x)**2) + (($y - $patch_y)**2));
if ($distance <= $max_distance) {
$nearby_patches{"${distance} $_"} = $_;
}
if ($verbose > 3) {
print " [", $_, "] --> ", $distance, "\n";
}
}
my $patch_count = 0;
my @stat_sum = (0,0,0,0,0,0);
my $weight_sum = 0;
foreach (sort {$a <=> $b} keys %nearby_patches) {
$patch_count = $patch_count + 1;
my $distance = $_;
my $stat_string = $all_patches{$nearby_patches{$_}};
if ($verbose > 2) {
print " ", $distance, ": ", $nearby_patches{$_}, " = ", $stat_string, "\n";
}
if ($distance == 0) {
$pixel_stats{"${x} ${y}"} = $stat_string;
if ($verbose > 2) {
print " gets ", $stat_string, "\n";
} else {
if ($verbose > 1) {
print "[", $x, ", ", $y, "] -> ", $pixel_stats{"${x} ${y}"}, "\n";
}
}
next PIXEL;
}
my $weight = (1.0 / ($distance));
my @stats = split " ", $stat_string;
for (my $stat = 0; $stat < 6; $stat++) {
$stat_sum [$stat] = $stat_sum [$stat] + ($stats [$stat] * $weight);
}
$weight_sum = $weight_sum + $weight;
}
for (my $stat = 0; $stat < 6; $stat++) {
$stat_sum [$stat] = $stat_sum [$stat] / $weight_sum;
}
$pixel_stats {"${x} ${y}"} = join " ", @stat_sum;
if ($verbose > 2) {
print " gets ", $pixel_stats{"${x} ${y}"}, "\n";
} else {
if ($verbose > 1) {
print "[", $x, ", ", $y, "] -> ", $pixel_stats{"${x} ${y}"}, "\n";
}
}
}
}
if ($verbose > 0) {
print "Generating a ", $image_width, " by ", $image_height, " image.\n";
}
for (my $y = 0; $y < $image_height; $y++) {
for (my $x = 0; $x < $image_width; $x++) {
my ($red_mean, $green_mean, $blue_mean, $red_stddev, $green_stddev,
$blue_stddev) = split " ", $pixel_stats{"${x} ${y}"};
if ($verbose > 3) {
print "[", $x, ", ", $y, "] -> ", $red_mean, ", ", $green_mean, ", ", $blue_mean, "\n";
}
my $red_val = &random_component_value ($red_mean, $red_stddev);
my $green_val = &random_component_value ($green_mean, $green_stddev);
my $blue_val = &random_component_value ($blue_mean, $blue_stddev);
$red_val = &to_srgb ($red_val);
$green_val = &to_srgb ($green_val);
$blue_val = &to_srgb ($blue_val);
$red_val = &roundoff ($red_val);
$green_val = &roundoff ($green_val);
$blue_val = &roundoff ($blue_val);
my $pixel_color = $out_image->colorAllocate ($red_val, $green_val, $blue_val);
$out_image->setPixel($x, $y, $pixel_color);
}
}
} else {
if ($verbose > 0) {
print "Generating a ", $image_width, " by ", $image_height, " image.\n";
}
while (my ($position, $stat_string) = each %all_patches) {
my ($left_x, $top_y) = split " ", $position;
my ($red_mean, $green_mean, $blue_mean, $red_stddev, $green_stddev,
$blue_stddev) = split " ", $stat_string;
if ($verbose > 3) {
print "[", $left_x, ", ", $top_y, "]:[", $left_x + $patch_size - 1, ", ",
$top_y + $patch_size - 1, "] -> ", $red_mean, ", ", $green_mean, ", ",
$blue_mean, "\n";
}
for (my $x = $left_x; $x < ($left_x + $patch_size); $x++) {
for (my $y = $top_y; $y < ($top_y + $patch_size); $y++) {
my $red_val = &random_component_value ($red_mean, $red_stddev);
my $green_val = &random_component_value ($green_mean, $green_stddev);
my $blue_val = &random_component_value ($blue_mean, $blue_stddev);
$red_val = &to_srgb ($red_val);
$green_val = &to_srgb ($green_val);
$blue_val = &to_srgb ($blue_val);
$red_val = &roundoff ($red_val);
$green_val = &roundoff ($green_val);
$blue_val = &roundoff ($blue_val);
my $pixel_color = $out_image->colorAllocate ($red_val, $green_val, $blue_val);
$out_image->setPixel($x, $y, $pixel_color);
}
}
}
}
#
# Write the file.
#
open OUTFILE, ">", $out_name;
binmode (OUTFILE);
print OUTFILE $out_image->png;
close OUTFILE;
#
# Given a mean and standard deviation, generate a component value.
# The value is clamped to (0,8192) to keep it valid in scRGB.
#
sub random_component_value {
my ($mean, $std_dev) = @_;
my $value = random_normal (1, $mean, $std_dev);
$value = 0 if ($value < 0);
$value = 8192 if ($value > 8192);
return $value;
}
#
# Round off a value to the nearest integer.
#
sub roundoff {
my ($value) = @_;
my $rounded_value = int($value + 0.5);
return $rounded_value;
}
#
# Convert an sRGB value to scRGB, based on IEC 61966-2-2.
#
sub to_scrgb {
my ($value) = @_;
my $result;
if ($value < 10.0) {
$result = $value * 2.4865;
return $result;
}
$result = ($value + 14.025) / 269.025;
$result = $result**2.4;
$result = $result * 8192.0;
return $result;
}
#
# Convert an scRGB value to sRGB, based on IEC 61966-2-2.
#
sub to_srgb {
my ($value) = @_;
my $result;
if ($value < 0.0) {
$result = 0.0;
return $result;
}
$result = $value / 8192.0;
if ($result < 0.00304) {
$result = 255.0 * (12.92 * $result);
return $result;
}
if ($result <= 1.0) {
$result = $result ** (1.0 / 2.4);
$result = $result * 269.025;
$result = $result - 14.025;
return $result;
}
$result = 255.0;
return $result;
}
__END__
=head1 NAME
fill_erasures - fill in erased parts of an image
=head1 SYNOPSIS
fill_erasures [options] input_file output_file
The input and output files are pictures in PNG format.
Options:
--help brief help message
--man full documentation
--verbose output progress information, can be repeated
--patch_size size of patches to process
--blend_patches if specified, blend adjacent patches together
=head1 OPTIONS
=over 8
=item B<help>
Print a brief help message and exit.
=item B<man>
Print the manual page and exit.
=item B<verbose>
Print informative progress and diagnostic messages to standard output.
The option can be repleated for additional output. High levels of verbosity
are intended for debugging the algorithms, and are somewhat cryptic.
=item B<patch_size>
The picture is devided into a number of patches, and the statistics for each
are computed to characterize their texture. The default for patch_size is 10,
meaning the patches are 10 by 10 pixels. Making the patch size larger will cause
the program to run faster, but the output won't be as smooth.
=item B<blend_patches>
Normally, each pixel has the statistical properties of its patch. If adjacent
patches are similar, this works well. If they are not, you can increase the
smoothness of the output by specifying the blend_patches option, which causes
each pixel's statistics to be icomputed from its own and the nearby patches.
Doing patch blending takes considerable time on large pictures.
=back
=head1 DESCRIPTION
Suppose you have a
picture which includes objects you wish to erase, revealing more of the sky
or other textured background. Unless the texture is very simple you cannot
copy it from elsewhere in the picture, since it will not match. You need
a realistic texture to replace the erased objects. fill_erasures will create
such a texture. Here is how to use it.
Erase any objects you wish to remove from the picture. Make sure the erased
areas are surrounded by the texture you wish to substitute for the objects.
Have fill_erasures read your picture as its input. The output file will have
erasures replaced by their surrounding textures. Use this output file as a layer
beneath your picture, so the replacement texture will show in place of the
erased areas.
=head1 BACKGROUND
The color of the sky is remarkably complex. Each small patch is not a single color, but
a small volume of colors. That distribution needs to be preserved or the synthesized sky
looks like a cartoon sky. In addition, the mean color changes with the angle from the
horizon. All three color components increase in intensity quickly as you begin looking
upward from
the horizon, then slowly decrease. Blue at first levels off at its maximum intensity,
then slowly declines. Green begins falling more quickly than blue, and falls at a steeper
angle. Red begins falling more quickly than green, and likewise falls faster at first,
but then begins to level off at high angles. If the photograph includes both tall and
short objects at a distance, this horizon effect can be observed to follow the contour
of the effective horizon.
=cut