def _extract(self):
        c_a = numpy.zeros((numpy.divide(self.container.shape[0], 2),
                           numpy.divide(self.container.shape[1], 2), self.container.shape[2]),
                          dtype=numpy.single)
        c_h = numpy.copy(c_a)
        c_v = numpy.copy(c_a)
        c_d = numpy.copy(c_a)

        self.watermark = numpy.zeros(self.aux[EmbeddingMethodConfig.AUX_STEGO_SIZE], dtype=self.container.dtype)

        for i in range(self.container.ndim):
            c_a[:, :, i], (c_h[:, :, i], c_v[:, :, i], c_d[:, :, i]) = pywt.dwt2(self.container[:, :, i],
                                                                                 self.aux[DeyMethodConfig.AUX_WAVELET])

            caw, (chw, cvw, cdw) = pywt.dwt2(self.stego[:, :, i], self.aux[DeyMethodConfig.AUX_WAVELET])

            cas = (caw - (1 - self.aux[DeyMethodConfig.AUX_G]) * c_a[:, :, i]) / self.aux[DeyMethodConfig.AUX_G]
            chs = (chw - (1 - self.aux[DeyMethodConfig.AUX_G]) * c_h[:, :, i]) / self.aux[DeyMethodConfig.AUX_G]
            cvs = (cvw - (1 - self.aux[DeyMethodConfig.AUX_G]) * c_v[:, :, i]) / self.aux[DeyMethodConfig.AUX_G]
            cds = (cdw - (1 - self.aux[DeyMethodConfig.AUX_G]) * c_d[:, :, i]) / self.aux[DeyMethodConfig.AUX_G]

            size = cas.shape[0] / 2

            LL = (cas[0:size, 0:size] + chs[0:size, 0:size] + cvs[0:size, 0:size] + cds[0:size, 0:size]) / 4
            LH = (cas[size:2 * size, 0:size] + chs[size:2 * size, 0:size] + cvs[size:2 * size, 0:size] + cds[
                                                                                                         size:2 * size,
                                                                                                         0:size]) / 4
            HL = (cas[0:size, size:2 * size] + chs[0:size, size:2 * size] +
                  cvs[0:size, size:2 * size] + cds[0:size, size:2 * size]) / 4
            HH = (cas[size:2 * size, size:2 * size] + chs[size:2 * size, size:2 * size] +
                  cvs[size:2 * size, size:2 * size] + cds[size:2 * size, size:2 * size]) / 4

            self.watermark[:, :, i] = pywt.idwt2((LL, (LH, HL, HH)), self.aux[DeyMethodConfig.AUX_WAVELET])
        super(DeyExtractor, self)._extract()
        return self.watermark

def collect(S, wavelet, mode, level):
    '''
    Returns the full quad tree of wavelet packets.
    @param S:         Input signal.
                      Both single and double precision floating-point data types are supported
                      and the output type depends on the input type. If the input data is not
                      in one of these types it will be converted to the default double precision
                      data format before performing computations.
    @param wavelet:   Wavelet to use in the transform. 
                      This must be a name of the wavelet from the wavelist() list.
    @param mode:      Signal extension mode to deal with the border distortion problem.
    @param level:     Number of decomposition steps to perform. If the level is None, then the
                      full decomposition up to the level computed with dwt_max_level() function for
                      the given data and wavelet lengths is performed.
    @return:          The full quad tree of wavelet packets.
    '''
    Nodes = [[] for i in range(level)]
    (CA, (CH, CV, CD)) = pywt.dwt2(S, wavelet=wavelet, mode=mode)
    Nodes[0] = [node.Node(CA, 0, 0), node.Node(CH, 0, 1), node.Node(CV, 0, 2), node.Node(CD, 0, 3)]
    for l in range(0, level-1):
        Parents = Nodes[l]
        Childs = []
        for p in range(len(Parents)):
            (CA, (CH, CV, CD)) = pywt.dwt2(Parents[p].C, wavelet=wavelet, mode=mode)
            Childs.append(node.Node(CA, l+1, 4*p))
            Childs.append(node.Node(CH, l+1, 4*p+1))
            Childs.append(node.Node(CV, l+1, 4*p+2))
            Childs.append(node.Node(CD, l+1, 4*p+3))
        Nodes[l+1] = Childs 
    return Nodes

def image_diff_dwt(lhs_image, rhs_image) -> int:
    _, (lhs_LH, lhs_HL, lhs_HH) = pywt.dwt2(lhs_image, 'haar')
    _, (rhs_LH, rhs_HL, rhs_HH) = pywt.dwt2(rhs_image, 'haar')
    d1 = cv2.absdiff(lhs_LH, rhs_LH).sum()
    d2 = cv2.absdiff(lhs_HL, rhs_HL).sum()
    d3 = cv2.absdiff(lhs_HH, rhs_HH).sum()
    return d1 + d2 + d3

def waveletTransformFunction(img1, img2):
    # whole of the wavelet families: ['haar', 'db', 'sym', 'coif', 'bior', 'rbio', 'dmey']
    imgVl = cv2.cvtColor(img1, cv2.COLOR_RGB2GRAY)
    wImgV= pywt.dwt2(imgVl,'haar')
    cAv, (cHv, cVv, cDv) = wImgV
    print (imgVl.shape)
    imgIr = cv2.cvtColor(img2, cv2.COLOR_RGB2GRAY)
    wImgIr = pywt.dwt2(imgNir, 'haar')
    cAi, (cHi, cVi, cDi) = wImgIr
    # cv2.imwrite('temporal.png', cAv)
    # cv2.imshow('cAv', cAv)
    # cv2.imshow('cHv', cHv)
    # cv2.imshow('cVv', cVv)
    # cv2.imshow('cDv', cDv)
    # cv2.imshow('imgVl', imgVl)
    # cAv1 = np.int16(cAv)
    # cHv1 = np.int16(cHv)
    # cVv1 = np.int16(cVv)
    # cDv1 = np.int16(cDv)
    # np.savetxt('texcAv.txt', cAv1)
    # np.savetxt('texcHv.txt', cHv)
    # print (cHv1[16:20][:18])
    # a short code for pyramid
    imgVlPy = []
    imgVlPy = cv2.pyrDown(imgVl,)
    #print (imgVlPy[16:20][:18])
    cv2.imshow('pyramid image', imgVlPy)
    cv2.waitKey(0)
    cv2.destroyAllWindows()
    a = imgVlPy[0:19][0:17]
    print (a)
    print a.shape

    def _embed(self):
        container = EmbeddingMethodStack.rgb_2_ycbcr(self.container)
        stego = numpy.copy(container)
        self.watermark = skimage.color.rgb2gray(self.watermark)
        self.watermark = scipy.misc.imresize(self.watermark, (numpy.divide(container.shape[0], 64),
                                                              numpy.divide(container.shape[1], 64)), interp='bicubic')
        super(ElahianEmbedder, self)._embed()

        watermark_bitstream = EmbeddingMethodStack.matrix_2_bitstream(self.watermark.astype(numpy.uint8))

        position = EmbeddingMethodStack.pseudo_rand_mask_create(numpy.asarray((numpy.divide(container.shape[0], 8),
                                                                               numpy.divide(container.shape[1], 8))))
        self.aux[ElahianMethodConfig.AUX_POSITION] = position

        c1a, (c1h, c1v, c1d) = pywt.dwt2(container[:, :, 0], self.aux[ElahianMethodConfig.AUX_WAVELET])
        c2a, (c2h, c2v, c2d) = pywt.dwt2(c1a, self.aux[ElahianMethodConfig.AUX_WAVELET])
        c3a, (c3h, c3v, c3d) = pywt.dwt2(c2a, self.aux[ElahianMethodConfig.AUX_WAVELET])

        x = position[:, 1]
        y = position[:, 0]
        c3ae = numpy.copy(c3a)
        # print c3ae.shape, watermark_bitstream.shape, position.shape
        for k in range(watermark_bitstream.size - 1):
            c3ae[y[k], x[k]] = c3a[y[k], x[k]] + self.aux[ElahianMethodConfig.AUX_G] * watermark_bitstream[k]

        c2ae = pywt.idwt2((c3ae, (c3h, c3v, c3d)), self.aux[ElahianMethodConfig.AUX_WAVELET])
        c1ae = pywt.idwt2((c2ae, (c2h, c2v, c2d)), self.aux[ElahianMethodConfig.AUX_WAVELET])
        stego[:, :, 0] = pywt.idwt2((c1ae, (c1h, c1v, c1d)), self.aux[ElahianMethodConfig.AUX_WAVELET])

        stego = EmbeddingMethodStack.ycbcr_2_rgb(stego, self.container.dtype)

        # print stego.shape, self.watermark.shape, self.container.shape

        return stego

 def getGeneralStatistics(self, hara=False, zern=False, tamura=False, only1D=None):
     generalStatistics = []
     if self.rows == 1 and self.columns == 1:
         for index in range(3):
             generalStatistics.append(self.image[0, 0, index])
         return generalStatistics
     if not only1D is None:
         im = only1D
         generalStatistics.extend(self._calculateStatistics(im, haralick=hara, zernike=zern))
         fourierTransform = np.abs(fftpack.fft2(im))  # fourierTransform
         generalStatistics.extend(self._calculateStatistics(fourierTransform))
         waveletTransform = pywt.dwt2(im, "sym5")[0]
         generalStatistics.extend(self._calculateStatistics(waveletTransform))
         waveletFourierTransform = pywt.dwt2(fourierTransform, "sym5")[0]
         generalStatistics.extend(self._calculateStatistics(waveletFourierTransform))
         if tamura:
             generalStatistics.extend(self.get3Dstatistics(tamura=True))
         return generalStatistics
     for index in range(3):
         im = self.image[:, :, index]
         generalStatistics.extend(self._calculateStatistics(im, haralick=hara, zernike=zern))
         fourierTransform = np.abs(fftpack.fft2(im))  # fourierTransform
         generalStatistics.extend(self._calculateStatistics(fourierTransform))
         waveletTransform = pywt.dwt2(im, "sym5")[0]
         generalStatistics.extend(self._calculateStatistics(waveletTransform))
         waveletFourierTransform = pywt.dwt2(fourierTransform, "sym5")[0]
         generalStatistics.extend(self._calculateStatistics(waveletFourierTransform))
     if tamura:
         generalStatistics.extend(self.get3Dstatistics(tamura=True))
     return generalStatistics

    def _extract(self, img):
        cA, (cH, cV, cD) = dwt2(img, self.mother)
        cH2 = cH.reshape(cH.size)
        cV2 = cV.reshape(cV.size)
        cD2 = cD.reshape(cD.size)
        assert cH2.shape == cV2.shape == cD2.shape
        buffers = (cH2, cV2, cD2)
        chunk_size = (cH2.size * 3) // (self.total_bits)
        seq0 = self.seq0[:chunk_size]
        seq1 = self.seq1[:chunk_size]

        byte = 0
        output = bytearray()
        for i in range(self.total_bits):
            target = buffers[i % 3]
            offset = (i//3) * chunk_size
            chunk = target[offset : offset + seq0.size]
            corr0, _ = pearsonr(chunk, seq0)
            corr1, _ = pearsonr(chunk, seq1)
            bit = int(corr1 > corr0)
            byte = (byte << 1) | bit
            if i % 8 == 7:
                output.append(byte)
                byte = 0
        #print repr(output)
        return output

def wavelet_fusion(Im1, Im2):
    m1 = np.mean(Im1)
    m2 = np.mean(Im2)
    s1 = np.std(Im1)
    s2 = np.std(Im2)
    g = s2/s1
    offset = m2-g*m1
    ImH = Im1*g+offset
    # Wavelet Transform step
    coeffs1 = pywt.dwt2(ImH, 'db1')
    coeffs2 = pywt.dwt2(Im2, 'db1')
    # SELECTION OF COEFFICIENT USING Linear combination BETWEEN both approximations.
    A = 2*(coeffs1[0]+coeffs2[0])/2
    # INVERSE WAVELET TRANSFORM TO GENERATE THE FUSED IMAGE.
    Xsyn = pywt.idwt2((A,coeffs1[1]), 'db1')
    return Xsyn

def recompose(image, coarse, wtype):

   print image.shape
   
   # nothing to be done: case in which coarse is the size of image
   if coarse.shape == image.shape:
      return coarse

   coeffs2 = pywt.dwt2(image,wtype,axes=(0,1))
   if coeffs2[0].shape == coarse.shape:
      ncoarse = coarse*2.0
   elif coeffs2[0].shape[0] > coarse.shape[0] and coeffs2[0].shape[1] > coarse.shape[1]:
      ncoarse = recompose(coeffs2[0],coarse*2.0,wtype)
   else:
      print "ERROR: Are you sure that \'%s\' is the right wavelet? Sizes don't match!"%wtype
      print image.shape
      print coarse.shape
      exit()
      

   #adjust size of upsampled version
   isz = coeffs2[0].shape
   ncoarse = ncoarse[0:isz[0],0:isz[1],:]

   ncoeffs2 = (ncoarse, coeffs2[1])
   ret = pywt.idwt2(ncoeffs2,wtype,axes=(0,1))

   return ret

 def applyTwoDDWTUsingLibrary(self, blockList):
     '''2D DWT using the pywt library.'''
     coeffs = pywt.dwt2(blockRef.getValue(), 'haar')
     cA, (cH, cV, cD) = coeffs
     print cA
     print cV
     print coeffs

def recompose(prefix, levels, suffix, wtype, cur=0):

   print prefix+str(cur)+suffix
   image = piio.read(prefix+str(cur)+suffix)
   print image.shape

   if levels==1: # if last level
      return image

   LL = 2 * recompose(prefix,levels-1,suffix, wtype, cur+1) ## 2 compensates the factor used in decompose 

   coeffs2  = pywt.dwt2(image,wtype,axes=(0,1))

   if coeffs2[0].shape != LL.shape:
      print "ERROR: Are you sure that \'%s\' is the right wavelet? Sizes don't match!"%wtype
      exit()

   ncoeffs2 = (LL, coeffs2[1])

   ret = pywt.idwt2(ncoeffs2,wtype, axes=(0,1))

   # adjust size of the upsampled image
   ret = ret[:image.shape[0],:image.shape[1],:] 
   print ret.shape
   return ret

    def _embed(self, img, payload, k):
        cA, (cH, cV, cD) = dwt2(img, self.mother)
        #assert cH2.shape == cV2.shape == cD2.shape
        buffer = numpy.zeros(cH.size + cV.size + cD.size)
        i = 0
        for arr in (cH, cV, cD):
            for row in arr:
                buffer[i:i+row.size] = row
                i += row.size
        chunk_size = buffer.size // self.total_bits
        sequence_of = (self.seq0[:chunk_size], self.seq1[:chunk_size])
        
        for i, bit in enumerate(iterbits(payload)):
            seq = sequence_of[bit]
            buffer[i * chunk_size : i * chunk_size + seq.size] += k * seq
        
        detail = (numpy.zeros(cH.shape), numpy.zeros(cV.shape), numpy.zeros(cD.shape))
        i = 0
        for arr in detail:
            h, w = arr.shape
            for r in range(h):
                arr[r] = buffer[i:i+w]
                i += w

        h, w = img.shape
        return idwt2((cA, detail), self.mother)[:h,:w]

    def _extract(self, img):
        cA, (cH, cV, cD) = dwt2(img, self.mother)
        cH2 = cH.reshape(cH.size)
        cV2 = cV.reshape(cV.size)
        assert cH2.shape == cV2.shape
        chunk_size = cH2.size // (self.total_bits // 2)
        seq0 = self.seq0[:chunk_size]
        seq1 = self.seq1[:chunk_size]

        byte = 0
        output = bytearray()
        for i in range(self.total_bits):
            target = (cH2, cV2)[i % 2]
            offset = (i//2) * chunk_size
            chunk = target[offset : offset + seq0.size]
            #chunk = target[i::self.total_bits][:seq0.size]
            #if not all(chunk[i] == chunk[i+1] for i in range(chunk.size-1)):
            corr0, _ = pearsonr(chunk, seq0)
            corr1, _ = pearsonr(chunk, seq1)
            bit = int(corr1 > corr0)
            #else:
            #    bit = 0
            byte = (byte << 1) | bit
            if i % 8 == 7:
                output.append(byte)
                byte = 0
        print repr(output)
        return output

    def _decompose(self, img):
        """
        Decompose an image into the WSQ subband pattern.
        
        Parameters
        ----------
        img : numpy array holding the image to be decomposed
        
        Returns
        -------
        subbands : list of 64 numpy arrays containing the WSQ subbands in order 
        """
        wavelet='coif1'
        subbands = []
        # first decompose image into 16 subbands
        temp1 = self._decompose16(img, wavelet)

        # next, decompose top left three subbands again into 16
        temp2 = []
        for i in xrange(3):
            temp2.append(self._decompose16(temp1[i], wavelet))

        # finally, decompose top left subband again into 4
        ll, hvd = pywt.dwt2(temp2[0][0], wavelet, mode='per')

        # insert subbands into list in correct order
        subbands.append(ll)
        subbands.extend(hvd)
        subbands.extend(temp2[0][1:])
        subbands.extend(temp2[1])
        subbands.extend(temp2[2])
        subbands.extend(temp1[3:])
        return subbands

def WaveletTransform(image, wavelet):
    coeffs = pywt.dwt2(image, wavelet)
    cA, (cH, cV, cD) = coeffs
    
    #####------For db2 trimming down the image to maintain size consistancy------#####
    if wavelet=='db2':
        cA=cA[1:len(cA),1:len(cA)]
        cH=cH[1:len(cH),1:len(cH)]
        cV=cV[1:len(cV),1:len(cV)]
        cD=cD[1:len(cD),1:len(cD)]
    print 'len(cA) '+str(len(cA))
    print 'len(cH) '+str(len(cH))
    print 'len(cV) '+str(len(cV))
    print 'len(cD) '+str(len(cD))
    
    #####------Scaling the transformed image by 2------#####
    cA=cv2.pyrUp(cA)
    cH=cv2.pyrUp(cH)
    cV=cv2.pyrUp(cV)
    cD=cv2.pyrUp(cD)
   
    print 'len(cA)up '+str(len(cA))
    print 'len(cH)up '+str(len(cH))
    print 'len(cV)up '+str(len(cV))
    print 'len(cD)up '+str(len(cD))
   
    return cA,cH,cV,cD

def test_per_axis_wavelets_and_modes():
    # tests seperate wavelet and edge mode for each axis.
    rstate = np.random.RandomState(1234)
    data = rstate.randn(16, 16, 16)

    # wavelet can be a string or wavelet object
    wavelets = (pywt.Wavelet('haar'), 'sym2', 'db4')

    # mode can be a string or a Modes enum
    modes = ('symmetric', 'periodization',
             pywt._extensions._pywt.Modes.reflect)

    coefs = pywt.dwtn(data, wavelets, modes)
    assert_allclose(pywt.idwtn(coefs, wavelets, modes), data, atol=1e-14)

    coefs = pywt.dwtn(data, wavelets[:1], modes)
    assert_allclose(pywt.idwtn(coefs, wavelets[:1], modes), data, atol=1e-14)

    coefs = pywt.dwtn(data, wavelets, modes[:1])
    assert_allclose(pywt.idwtn(coefs, wavelets, modes[:1]), data, atol=1e-14)

    # length of wavelets or modes doesn't match the length of axes
    assert_raises(ValueError, pywt.dwtn, data, wavelets[:2])
    assert_raises(ValueError, pywt.dwtn, data, wavelets, mode=modes[:2])
    assert_raises(ValueError, pywt.idwtn, coefs, wavelets[:2])
    assert_raises(ValueError, pywt.idwtn, coefs, wavelets, mode=modes[:2])

    # dwt2/idwt2 also support per-axis wavelets/modes
    data2 = data[..., 0]
    coefs2 = pywt.dwt2(data2, wavelets[:2], modes[:2])
    assert_allclose(pywt.idwt2(coefs2, wavelets[:2], modes[:2]), data2,
                    atol=1e-14)

def htrans(A):
    h0 = A
    res = []
    while h0.shape[0]>1 and h0.shape[1]>1:
        h0, (hx, hy, hc) = pywt.dwt2(h0, 'haar')
        res = [(h0, h0, h0)]+res
    out, _ = pywt.coeffs_to_array([h0]+res, padding=1)
    return out

def test_idwt2_axes():
    data = np.array([[0, 1, 2, 3],
                     [1, 1, 1, 1],
                     [1, 4, 2, 8]])
    coefs = pywt.dwt2(data, 'haar', axes=(1, 1))
    assert_allclose(pywt.idwt2(coefs, 'haar', axes=(1, 1)), data, atol=1e-14)

    # too many axes
    assert_raises(ValueError, pywt.idwt2, coefs, 'haar', axes=(0, 1, 1))

 def _embed(self, img, payload, k):
     rand = Random(self.seed)
     cA, (cH, cV, cD) = dwt2(img, self.mother)
     cD2 = cD.reshape(cD.size)
     for bit in iterbits(payload):
         seq = numpy.array([int(rand.random() > 0.95) for _ in range(cD2.size)]) 
         if bit:
             cD2 += k * seq
     return idwt2((cA, (cH, cV, cD2.reshape(cD.shape))), self.mother)[:img.shape[0],:img.shape[1]]

def test_ignore_invalid_keys():
    data = np.array([[0, 4, 1, 5, 1, 4], [0, 5, 6, 3, 2, 1], [2, 5, 19, 4, 19, 1]])

    wavelet = pywt.Wavelet("haar")

    LL, (HL, LH, HH) = pywt.dwt2(data, wavelet)
    d = {"aa": LL, "da": HL, "ad": LH, "dd": HH, "foo": LH, "a": HH}

    assert_allclose(pywt.idwt2((LL, (HL, LH, HH)), wavelet), pywt.idwtn(d, wavelet), atol=1e-15)