Let's say I have list of genomic positions, e.g. SNP calls in a VCF file, and a BAM file with read alignments. What would be the best way to, for each of those positions, determine if in the BAM file there are reads with indels either at/overlapping the position, or in the vicinity of the position, say 5 bp upstream or downstream. This is a read level query, it's not indel calling. Basically I want to for a given position consider the ±5 bp window around the position and ask if there are any indels within this window, at the level of reads, and if so count the number of indels in the window.
I think the CIGAR strings of the BAM file probably holds all the information needed, but the question is how to do this in a simple, practical way.