What you're looking to do is fairly straightforward. I'd just call SNPs with Samtools, then merge the lists with a little script.

From that look of it, what you are trying to do is the standard SNP calling procedure. You might want to try out the bwa -> GATK pipeline where you can find: https://www.broadinstitute.org/gatk/

Once you've got the vcf, it will be easy to obtain the information you want (something like homozygous variation on the sample)