Calculate a square matrix of r2 between all SNPs in genome during LD analysis

Question

[Denis]

I have a multiple genotypes in vcf format. I'd like to calculate a square matrix of r2 between each pair of SNPs in genome (even for these located in differen chromosomes). Below is a toy example i wish to get:

       rs1   rs2   rs3
rs1   1      0.8   0.3
rs2   0.8    1     0.5
rs3   0.3    0.5   1

Besides, how can i handle observations in which distance between two SNPs will be NA ?

Many thanks, Denis