shadow_mask will preserve previous alpha

This makes sure that one can apply shadow_mask after grassfirealpha. (The opposite is also possible but has a different effect).

src/asp/PhotometryTK/phoitalbedo.cc

@@ -253,6 +253,8 @@ int main( int argc, char *argv[] ) {
     // Setting up working level
     if ( opt.level < 0 )
       opt.level = drg_plate->num_levels() - 1;
+    if ( opt.level >= drg_plate->num_levels() )
+      vw_throw( ArgumentErr() << "Can't request level higher than available in DRG plate" );
 
     // Diving up jobs and deciding work units
     std::list<BBox2i> workunits;
@@ -279,7 +281,8 @@ int main( int argc, char *argv[] ) {
     }
 
     // Double check for an iteration error
-    if ( albedo_plate->num_levels() == 0 ) {
+    if ( albedo_plate->num_levels() == 0 ||
+         opt.level >= albedo_plate->num_levels() ) {
       project_info.set_current_iteration(0);
       remote_ptk->UpdateIteration(0);
     }

src/asp/PhotometryTK/shadow_mask.cc

@@ -37,13 +37,28 @@ void shadow_mask_nodata( Options& opt,
                                          georef, TerminalProgressCallback("photometrytk","Writing:"));
 }
 
+template <class PixelT>
+class ThresholdAlphaFunctor : public UnaryReturnSameType {
+  typedef typename vw::PixelChannelType<PixelT>::type channel_type;
+  double m_threshold;
+public:
+  ThresholdAlphaFunctor( double t ) : m_threshold(t) {}
+
+  template <class ArgT>
+  inline ArgT operator()( ArgT const& arg ) const {
+    ArgT copy = arg;
+    if ( PixelGray<channel_type>(copy) <= m_threshold ) {
+      copy.a() = 0;
+    }
+    return copy;
+  }
+};
+
 template <class PixelT>
 void shadow_mask_alpha( Options& opt,
                         std::string input,
                         std::string output ) {
   typedef typename PixelChannelType<PixelT>::type ChannelT;
-  typedef typename PixelWithoutAlpha<PixelT>::type PNAlphaT;
-  typedef typename MaskedPixelType<PNAlphaT>::type PMaskT;
   if ( opt.threshold < 0 ) {
     ChannelRange<ChannelT> helper;
     opt.threshold = 0.1*helper.max();
@@ -52,10 +67,9 @@ void shadow_mask_alpha( Options& opt,
   cartography::GeoReference georef;
   cartography::read_georeference(georef, input);
   DiskImageView<PixelT> input_image(input);
-  ImageViewRef<PMaskT> masked_input = alpha_to_mask(input_image);
-  ImageViewRef<PMaskT> result =
-    intersect_mask(masked_input,create_mask(threshold(apply_mask(masked_input),opt.threshold)));
-  cartography::write_georeferenced_image(output, mask_to_alpha(result), georef,
+  ImageViewRef<PixelT> result =
+    UnaryPerPixelView<DiskImageView<PixelT>,ThresholdAlphaFunctor<PixelT> >(input_image,ThresholdAlphaFunctor<PixelT>(opt.threshold));
+  cartography::write_georeferenced_image(output, result, georef,
                                          TerminalProgressCallback("photometrytk","Writing:"));
 }