If you have annotated genome, you can calculate exon coverage with either HTSeqlink text which has a command-line script that looks like:

htseq-count your.sam exons.gff  > feature_counts.txt

you can also use bedtools to do the same, see the example there that looks like this:

samtools view -bf 0x2 reads.bam | bamToBed -i stdin | coverageBed -a stdin -b exons.bed > exons.bed.proper.coverage

to calculate the coverage of aligned and paired reads over exons.

For merging, unless your aligner requires it, that's not necessary. You can do QC on the ends separately (and see large differences). For trimming and then discarding reads/pairs, you do want to make sure you discard from both or neither ends as nearly all aligners assume the fastq ands are in the same order.

Brent had a blog post on filtering paired end reads, could be a good start:

http://hackmap.blogspot.com/2010/09/filtering-paired-end-reads-high.html

As for the coverage question. One approach is to reconstitute the fragments from the read pairs and use these to compute the coverage.

filter paired end reads, you may try thishttps://github.com/Czh3/NGSTools/blob/master/qualityControl.py.

usage: qualityControl.py [-h] -1 FASTQ1 [-2 FASTQ2] [-q MINQUALITY]
                         [-p MINPERCENT] [-a AVERAGEQUALITY] [-f {33,64}] -o1
                         OUTFILE1 [-o2 OUTFILE2]

Illumina Single-end or Pair-end sequencing quality control
tools.<zhangchao3@hotmail.com>

optional arguments:
  -h, --help            show this help message and exit
  -1 FASTQ1, --fastq1 FASTQ1
                        The input fastq1 file.
  -2 FASTQ2, --fastq2 FASTQ2
                        The input fastq2 file.
  -q MINQUALITY, --minQuality MINQUALITY
                        Minimum quality score to keep.
  -p MINPERCENT, --minPercent MINPERCENT
                        Minimum percent of bases that must have [-q] quality.
  -a AVERAGEQUALITY, --averageQuality AVERAGEQUALITY
                        The average quality of SE/PE reads less than the
                        averageQuality will be descard.
  -f {33,64}, --phredQualtiy {33,64}
                        phred quality score
  -o1 OUTFILE1, --outFile1 OUTFILE1
                        fastq1 output file.
  -o2 OUTFILE2, --outFile2 OUTFILE2
                        fastq2 output file.

Heng Li has recently added new functionality to samtools (0.1.15) that will compute coverage based on a BED file ("-L" option). I believe it uses a BAM index file to seek to the appropriate regions via a disk seek. Depending on your disks, exons are probably near the point (N=~250,000 for knownGene transcripts) where it would be faster to use BEDTools coverageBed. My guess is that the speed plots for the two approaches intersect when the BED file has between 1e5 and 1e6 features.

Either way, samtools will use less memory if that is a concern. I don't think it tells you anything about the "shape" of the coverage, whereas coverageBed does.

This is a long way of telling you that you have lots of good options.

Qualimap computes coverage and many other QC issues from BAM file. Here's the link:

http://qualimap.bioinfo.cipf.es/