RuntimeWarning: Mean of empty slice ... or ... gdalbuildvrt ERROR 1: No input dataset specified.

[absudabsu]

Hi,

I am attempting to run s2p, but I am getting some runtime errors, and I am not sure how to proceed. This is what I've done: I acquired/assembled the XML and TIF files, and edited the config JSON with this information. There are only a few parameters I changed from the basic example that ships with the code: (testdata/input_pair/config.json).

1. I have three images, and no access to the original image that these were clipped from. So, if I am understanding correctly, I set 'full_img' = True
2. I updated the 'images' list with "img" and "rpc" for each image
3. I wasn't sure if 'utm_bbx' was used, but I populated this also with average "middle" values of the valid UTM ranges of each of the images, given the RPC parameters.

With this assumptions, I tried calling the s2p algorithm. It works great for the input_pair and input_triplet data, but when I run it on my own data (using 3-images from a dataset), I get the following runtime error:

WARNING: ignoring unknown parameter utm_bbx.                                                                                       
tile size: 286 286                                                                                                                 
total number of tiles: 98 (7 x 7) x 2 pairs                                                                                        
                                                                                                                                   
discarding masked tiles...                                                                                                         
done 49 / 49 tiles                                                                                                                 
Elapsed time: 0:00:00.136507                                                                                                       
                                                                                                                                   
correcting pointing locally...                                                                                                     
done 98 / 98 tiles                                                                                                                 
Elapsed time: 0:00:12.288325                                                                                                       
                                                                                                                                   
correcting pointing globally...                                                                                                    
Elapsed time: 0:00:00.016108                                                                                                       
                                                                                                                                   
rectifying tiles...                                                                                                                
Elapsed time: 0:00:03.823495                                                                                                       
                                                                                                                                   
running stereo matching...                                                                                                         
Elapsed time: 0:02:18.679599                                                                                                       
                                                                                                                                   
computing height maps...                                                                                                           
Elapsed time: 0:00:03.921784                                                                                                       
                                                                                                                                   
computing local pairwise height offsets...                                                                                         
done 49 / 49 tiles                                                                                                                 
Elapsed time: 0:00:00.218622                                                                                                       
                                                                                                                                   
computing global pairwise height offsets...                                                                                        
**./pipelines/1st/s2p.py:343: RuntimeWarning: Mean of empty slice**                                                                    
  global_mean_heights = np.nanmean(local_mean_heights, axis=0)                                                                     
merging height maps and computing point clouds...                                                                                  
Elapsed time: 0:00:02.440702                                                                                                       
                                                                                                                                   
computing global source window (xoff, yoff, xsize, ysize)...                                                                       
Elapsed time: 0:00:00.000659                                                                                                       
                                                                                                                                   
computing DSM by tile...                                                                                                           
Elapsed time: 0:00:00.114463                                                                                                       
                                                                                                                                   
computing global DSM...  

RUN: gdalbuildvrt -vrtnodata nan -input_file_list ./data/tmp/s2p_out/gdalbuildvrt_input_file_list.txt ./data/tmp/s2p_out/dsm.vrt   
ERROR 1: No input dataset specified.                                                                                               
Usage: gdalbuildvrt [-tileindex field_name]                                                                                        
                    [-resolution {highest|lowest|average|user}]                                                                    
                    [-te xmin ymin xmax ymax] [-tr xres yres] [-tap]                                                               
                    [-separate] [-b band] [-sd subdataset]                                                                         
                    [-allow_projection_difference] [-q]                                                                            
                    [-addalpha] [-hidenodata]                                                                                      
                    [-srcnodata "value [value...]"] [-vrtnodata "value [value...]"]                                                
                    [-a_srs srs_def]                                                                                               
                    [-r {nearest,bilinear,cubic,cubicspline,lanczos,average,mode}]                                                 
                    [-input_file_list my_list.txt] [-overwrite] output.vrt [gdalfile]*                                             
                                                                                                                                   
e.g.                                                                                                                               
  % gdalbuildvrt doq_index.vrt doq/*.tif                                                                                           
  % gdalbuildvrt -input_file_list my_list.txt doq_index.vrt                                                                        
                                                                                                                                   
NOTES:                                                                                                                             
  o With -separate, each files goes into a separate band in the VRT band.                                                          
    Otherwise, the files are considered as tiles of a larger mosaic.                                                               
  o -b option selects a band to add into vrt.  Multiple bands can be listed.                                                       
    By default all bands are queried.                                                                                              
  o The default tile index field is 'location' unless otherwise specified by                                                       
    -tileindex.                                                                                                                    
  o In case the resolution of all input files is not the same, the -resolution                                                     
    flag enable the user to control the way the output resolution is computed.                                                     
    Average is the default.                                                                                                        
  o Input files may be any valid GDAL dataset or a GDAL raster tile index.                                                         
  o For a GDAL raster tile index, all entries will be added to the VRT.                                                            
  o If one GDAL dataset is made of several subdatasets and has 0 raster bands,                                                     
    its datasets will be added to the VRT rather than the dataset itself.                                                          
    Single subdataset could be selected by its number using the -sd option.                                                        
  o By default, only datasets of same projection and band characteristics                                                          
    may be added to the VRT.                                                                                                       
Traceback (most recent call last):                                                                                                 
  File "./pipelines/1st/s2plib/common.py", line 90, in run                                                                         
    stderr=sys.stderr, env=env)                                                                                                    
  File "/opt/intel/intelpython35/lib/python3.5/subprocess.py", line 581, in check_call                                             
    raise CalledProcessError(retcode, cmd)                                                                                         
subprocess.CalledProcessError: Command 'gdalbuildvrt -vrtnodata nan -input_file_list ./data/tmp/s2p_out/gdalbuildvrt_input_file_lis
t.txt ./data/tmp/s2p_out/dsm.vrt' returned non-zero exit status 1                                                                  
                                                                                                                                   
During handling of the above exception, another exception occurred:   

Traceback (most recent call last):                                                                                                 
  File "<stdin>", line 1, in <module>                                                                                              
  File "<string>", line 225, in <module>                                                                                           
  File "<string>", line 194, in call_s2p_main                                                                                      
  File "./pipelines/1st/s2p.py", line 646, in main                                                                                 
    global_dsm(tiles)                                                                                                              
  File "./pipelines/1st/s2p.py", line 515, in global_dsm                                                                           
    out_dsm_vrt))                                                                                                                  
  File "./pipelines/1st/s2plib/common.py", line 97, in run                                                                         
    e.output})                                                                                                                     
s2plib.common.RunFailure: {'environment': environ({'TERM': 'xterm', 'JUPYTER_SERVER_ROOT': '/home/myuser', 'LINES': '62', 'PWD'
: '/home/myuser/projects/test-s2p', 'PATH': '/home/myuser/projects/test-s2p/pipelines/1st/bin:/opt/intel/intelpython35/
bin/:/opt/intel/intelpython35/bin/:/usr/local/sbin:/usr/local/bin:/usr/sbin:/usr/bin:/sbin:/bin', 'SHLVL': '1', '_': '/opt/intel/in
telpython35/bin/python', 'USER': 'myuser', 'PYXTERM_DIMENSIONS': '80x25', 'OLDPWD': '/home/myuser', 'JUPYTER_SERVER_URL': '
http://127.0.0.1:37281/user/myuser/', 'HOME': '/home/myuser', 'SHELL': '/bin/bash', 'COLUMNS': '131'}), 'output': None, 'co
mmand': 'gdalbuildvrt -vrtnodata nan -input_file_list ./data/tmp/s2p_out/gdalbuildvrt_input_file_list.txt ./data/tmp/s2p_out/dsm.vr
t'}

My JSON config_template looks like this:

{'clean_intermediate': False,                                                                                                      
 'clean_tmp': False,                                                                                                               
 'debug': False,                                                                                                                   
 'disable_srtm': True,                                                                                                             
 'disp_range_extra_margin': 0.2,                                                                                                   
 'disp_range_method': 'sift',                                                                                                      
 'disp_range_srtm_high_margin': 50,                                                                                                
 'disp_range_srtm_low_margin': -20,                                                                                                
 'dsm_resolution': 0.5,                                                                                                            
 'epipolar_thresh': 0.5,                                                                                                           
 'full_img': True,                                                                                                                 
 'fusion_thresh': 3,                                                                                                               
 'horizontal_margin': 20,                                                                                                          
 'images': [{'img': './data/tmp/s2p_in/s2p_INPUT_IM_0.tif',                                                                        
             'rpc': './data/tmp/s2p_in/s2p_INPUT_RPC_0.xml'},                                                                      
            {'img': './data/tmp/s2p_in/s2p_INPUT_IM_1.tif',                                                                        
             'rpc': './data/tmp/s2p_in/s2p_INPUT_RPC_1.xml'},                                                                      
            {'img': './data/tmp/s2p_in/s2p_INPUT_IM_2.tif',                                                                        
             'rpc': './data/tmp/s2p_in/s2p_INPUT_RPC_2.xml'}],                                                                     
 'matching_algorithm': 'mgm',                                                                                                      
 'msk_erosion': 0,                                                                                                                 
 'n_gcp_per_axis': 5,                                                                                                              
 'out_dir': './data/tmp/s2p_out/',                                                                                                 
 'roi': {},                                                                                                                        
 'sift_match_thresh': 0.6,                                                                                                         
 'skip_existing': False,                                                                                                           
 'subsampling_factor': 1,                                                                                                          
 'temporary_dir': './data/tmp/s2p_tmp/',                                                                                           
 'tile_size': 300,                                                                                                                 
 'utm_bbx': [352280.33713790023,                                                                                                   
             352580.33713790023,                                                                                                   
             6184003.9602460349,                                                                                                   
             6184303.9602460349],                                                                                                  
 'vertical_margin': 5} 

Any ideas? I was thinking that maybe the images don't cover the same area... which may cause the algorithm to be missing information for some locations. I'm not sure how to get around this, since it seems the utm_bbx was not even used, and it will be hard to tell before-hand whether or not two images contain information that s2p can exploit. Is my interpretation of the error correct? If yes: should it work if I call it with two "overlapping" images?

Thanks in advance!

[absudabsu]

I tried with a large batch of images, which are guaranteed to overlap. No luck, same error.

Looking at it again... Perhaps rather np.mean() returning NaN (a warning), the real error is: ERROR 1: No input dataset specified. ? I'm not sure how what I'm doing wrong.. Any hints?

[dyoussef]

Hello @sirgogo,
Thank you for bringing up the problem to our attention.
Please, can you send your dataset?
Thanks in advance

[absudabsu]

@dyoussef Thanks for the reply. The dataset I am using is available here: http://www.jhuapl.edu/satellite-benchmark.html . From my understanding, its the same dataset used in the IARPA top coder challenge, but I could be wrong. The specific images I am using are the ones available in the (above, specified) dataset, such as: iarpa_3D-satellite-dataset/Challenge_Data_and_Software/cropimagedata/Explorer/...

I did some more digging, and I think that maybe the issue is due to images not overlapping 100%, since I found that for some image pairs (just two at a time), the algorithm does run, but tries to overallocate memory (I see all kinds of malloc calls on the screen) when it is computing the global DSM. However, for other image pairs, it fails arbitrarily like above.

I tried a similar strategy on a smaller reconstruction are by setting 'full_img' = True and specifying a small ROI. For some image chips the reconstruction was able to complete (albeit on a very small area)! For other chips, it fails arbitrarily like above.

I have not inspected thoroughly, but from my initial experiments, it seems like the algorithm works/runs when two supplied images are very close (e.g. consecutive), such that their projections cover the same/similar geodesic bounding box with very high overlap, and their image-plane coordinates are also simliar.

[carlodef]

Hi sirgogo,
Could you please specify the filenames of the input images and provide your full json configuration file so that we can try to reproduce the failure?

From your output I see that the size of your input images is about 2000 x 2000 pixels (from 286 * 7). The original images of the IARPA challenge are much bigger. Maybe your issue comes from the fact that you've cropped the input images without modifying the corresponding RPC coefficients. There's a tool for that: rpc_cropper.