featureCounts can take multiple bam files as input: just use

featureCounts -parameters etc *.bam or $(ls *bam)

You'll get a counts file that is tab-delimited which you can then parse out

grep - v "#" counts_file | cut -d$'\t' -f1,7- > counts.matrix

If you just need to quantify read counts you can use a transcript aligner like Salmon instead, which goes straight from FASTQ files to count matrices (among a few other things), thus saving you the space of storing BAM files. Salmon tends to be quicker, more memory efficient, and more accurate than STAR + featureCounts anyway.

If you really want to use, or need to use featureCounts, you can save the individual count matrices, and merge/join all of them by gene name/id using your favorite programming language.

If you are running something separate on your bam, run RSEM. It's smarter.

But yes, Salmon or Kallisto would solve your space problem, and be faster.

You can try in bash script as well

# get the count
ls -1  *featureCount.txt | parallel 'cat {} | sed '1d' | cut -f7 {} > {/.}_clean.txt' 
ls -1  *featureCount.txt | head -1 | xargs cut -f1 > genes.txt
paste genes.txt *featureCount_clean.txt > output.txt