'Test that when auto-rescaling integer data with "blank" values (where
the blanks are replaced by NaN in the float data), that the "BLANK"
keyword is removed from the header.
Further, test that when using the ``scale_back=True`` option the blank
values are restored properly.
Regression test for https://github.com/astropy/astropy/issues/3865'

def test_scale_back_with_blanks(self):

    arr = np.arange(5, dtype=np.int32)
    hdu = fits.PrimaryHDU(data=arr)
    hdu.scale('int16', bscale=1.23)
    hdu.data[0] = 9999
    hdu.header['BLANK'] = 9999
    hdu.writeto(self.temp('test.fits'))
    with fits.open(self.temp('test.fits')) as hdul:
        data = hdul[0].data
        assert np.isnan(data[0])
        hdul.writeto(self.temp('test2.fits'))
    with catch_warnings() as w:
        with fits.open(self.temp('test2.fits')) as hdul2:
            assert (len(w) == 0)
            assert ('BLANK' not in hdul2[0].header)
            data = hdul2[0].data
            assert np.isnan(data[0])
    with fits.open(self.temp('test.fits'), scale_back=True, mode='update') as hdul3:
        data = hdul3[0].data
        assert np.isnan(data[0])
    with fits.open(self.temp('test.fits'), do_not_scale_image_data=True) as hdul4:
        assert (hdul4[0].header['BLANK'] == 9999)
        assert (hdul4[0].header['BSCALE'] == 1.23)
        assert (hdul4[0].data[0] == 9999)
