It's true these polynomials are not unimodal for n = 4, 5, 6. But they seem to be unimodal after that (I've just now calculated up to n=35). I had observed the non-unimodality of n=4,5,6 several years ago, but had thought it would continue to fail, so it's interesting to me that it seems to be unimodal at n=7 and beyond.
When the ASM is a permutation matrix, the right statistic is a weighted inversion count: weight each inversion pair by the larger number in the pair and then add up all the contributions. (I'm surprised I haven't added this permutation statistic to FindStat; I will remedy that!) The generating function for this statistic on permutations is an interesting analogue of the q-factorial. (See https://arxiv.org/pdf/1002.3391v2.pdf Corollary 6.)
I also think looking at the Gog or Magog triangles to try and find the right statistic is a good idea, but it's not one of the usual statistics, as far as I can tell. I've thought about this question quite a lot and would be very interested if anyone is able to make progress on finding this statistic on ASM or TSSCPP.