Python progressbar.ProgressBar類代碼示例

    def generate( self, out_path, aux, idx_in, idx_out ):
        scheme_high = amico.lut.create_high_resolution_scheme( self.scheme, b_scale = 1 )
        protocolHR = self.scheme2noddi( scheme_high )

        nATOMS = len(self.IC_ODs)*len(self.IC_VFs) + 1
        progress = ProgressBar( n=nATOMS, prefix="   ", erase=True )

        # Coupled contributions
        IC_KAPPAs = 1 / np.tan(self.IC_ODs*np.pi/2)
        for kappa in IC_KAPPAs:
            signal_ic = self.synth_meas_watson_SH_cyl_neuman_PGSE( np.array([self.dPar*1E-6, 0, kappa]), protocolHR['grad_dirs'], np.squeeze(protocolHR['gradient_strength']), np.squeeze(protocolHR['delta']), np.squeeze(protocolHR['smalldel']), np.array([0,0,1]), 0 )

            for v_ic in self.IC_VFs:
                dPerp = self.dPar*1E-6 * (1 - v_ic)
                signal_ec = self.synth_meas_watson_hindered_diffusion_PGSE( np.array([self.dPar*1E-6, dPerp, kappa]), protocolHR['grad_dirs'], np.squeeze(protocolHR['gradient_strength']), np.squeeze(protocolHR['delta']), np.squeeze(protocolHR['smalldel']), np.array([0,0,1]) )

                signal = v_ic*signal_ic + (1-v_ic)*signal_ec
                lm = amico.lut.rotate_kernel( signal, aux, idx_in, idx_out, False )
                np.save( pjoin( out_path, 'A_%03d.npy'%progress.i) , lm )
                progress.update()

        # Isotropic
        signal = self.synth_meas_iso_GPD( self.dIso*1E-6, protocolHR)
        lm = amico.lut.rotate_kernel( signal, aux, idx_in, idx_out, True )
        np.save( pjoin( out_path, 'A_%03d.npy'%progress.i) , lm )
        progress.update()

    def resample( self, in_path, idx_out, Ylm_out, doMergeB0 ) :
        if doMergeB0:
            nS = 1+self.scheme.dwi_count
            merge_idx = np.hstack((self.scheme.b0_idx[0],self.scheme.dwi_idx))
        else:
            nS = self.scheme.nS
            merge_idx = np.arange(nS)
        KERNELS = {}
        KERNELS['model'] = self.id
        KERNELS['D']     = np.zeros( (len(self.d_perps),181,181,nS), dtype=np.float32 )
        KERNELS['CSF']   = np.zeros( (len(self.d_isos),nS), dtype=np.float32 )

        nATOMS = len(self.d_perps) + len(self.d_isos)
        progress = ProgressBar( n=nATOMS, prefix="   ", erase=True )

        # Tensor compartment(s)
        for i in xrange(len(self.d_perps)) :
            lm = np.load( pjoin( in_path, 'A_%03d.npy'%progress.i ) )
            KERNELS['D'][i,...] = amico.lut.resample_kernel( lm, self.scheme.nS, idx_out, Ylm_out, False )[:,:,merge_idx]
            progress.update()

        # Isotropic compartment(s)
        for i in xrange(len(self.d_isos)) :
            lm = np.load( pjoin( in_path, 'A_%03d.npy'%progress.i ) )
            KERNELS['CSF'][i,...] = amico.lut.resample_kernel( lm, self.scheme.nS, idx_out, Ylm_out, True )[merge_idx]
            progress.update()

        return KERNELS

    def resample( self, in_path, idx_out, Ylm_out ):
        nATOMS = len(self.IC_ODs)*len(self.IC_VFs) + 1

        KERNELS = {}
        KERNELS['model'] = self.id
        KERNELS['wm']    = np.zeros( (nATOMS-1,181,181,self.scheme.nS), dtype=np.float32 )
        KERNELS['iso']   = np.zeros( self.scheme.nS, dtype=np.float32 )
        KERNELS['kappa'] = np.zeros( nATOMS-1, dtype=np.float32 )
        KERNELS['icvf']  = np.zeros( nATOMS-1, dtype=np.float32 )
        KERNELS['norms'] = np.zeros( (self.scheme.dwi_count, nATOMS-1) )

        progress = ProgressBar( n=nATOMS, prefix="   ", erase=True )

        # Coupled contributions
        for i in xrange( len(self.IC_ODs) ):
            for j in xrange( len(self.IC_VFs) ):
                lm = np.load( pjoin( in_path, 'A_%03d.npy'%progress.i ) )
                idx = progress.i - 1
                KERNELS['wm'][idx,:,:,:] = amico.lut.resample_kernel( lm, self.scheme.nS, idx_out, Ylm_out, False )
                KERNELS['kappa'][idx] = 1.0 / np.tan( self.IC_ODs[i]*np.pi/2.0 )
                KERNELS['icvf'][idx]  = self.IC_VFs[j]
                KERNELS['norms'][:,idx] = 1 / np.linalg.norm( KERNELS['wm'][idx,0,0,self.scheme.dwi_idx] ) # norm of coupled atoms (for l1 minimization)
                progress.update()

        # Isotropic
        lm = np.load( pjoin( in_path, 'A_%03d.npy'%progress.i ) )
        KERNELS['iso'] = amico.lut.resample_kernel( lm, self.scheme.nS, idx_out, Ylm_out, True )
        progress.update()

        return KERNELS

def debiasRician(DWI,SNR,mask,scheme):
    debiased_DWI = np.zeros(DWI.shape)
    t = time.time()
    progress = ProgressBar( n=mask.sum(), prefix="   ", erase=True )
    for ix in range(DWI.shape[0]):
        for iy in range(DWI.shape[1]):
            for iz in range(DWI.shape[2]):
                if mask[ix,iy,iz]:
                    b0 = DWI[ix,iy,iz,scheme.b0_idx].mean()
                    sigma_diff = b0/SNR
                    init_guess = DWI[ix,iy,iz,:].copy()
                    tmp = minimize(F_norm_Diff_K, init_guess, args=(init_guess,sigma_diff), method = 'L-BFGS-B', jac=der_Diff)
                    debiased_DWI[ix,iy,iz] = tmp.x
                    progress.update()
    print('   [ %s ]' % ( time.strftime("%Hh %Mm %Ss", time.gmtime(time.time()-t) ) ))
    return debiased_DWI

    def fit(self):
        """Fit the model to the data iterating over all voxels (in the mask) one after the other.
        Call the appropriate fit() method of the actual model used.
        """
        if self.niiDWI is None:
            raise RuntimeError('Data not loaded; call "load_data()" first.')
        if self.model is None:
            raise RuntimeError('Model not set; call "set_model()" first.')
        if self.KERNELS is None:
            raise RuntimeError(
                'Response functions not generated; call "generate_kernels()" and "load_kernels()" first.'
            )
        if self.KERNELS["model"] != self.model.id:
            raise RuntimeError("Response functions were not created with the same model.")

        self.set_config("fit_time", None)
        totVoxels = np.count_nonzero(self.niiMASK_img)
        print '\n-> Fitting "%s" model to %d voxels:' % (self.model.name, totVoxels)

        # setup fitting directions
        peaks_filename = self.get_config("peaks_filename")
        if peaks_filename is None:
            DIRs = np.zeros(
                [self.get_config("dim")[0], self.get_config("dim")[1], self.get_config("dim")[2], 3], dtype=np.float32
            )
            nDIR = 1
            gtab = gradient_table(self.scheme.b, self.scheme.raw[:, :3])
            DTI = dti.TensorModel(gtab)
        else:
            niiPEAKS = nibabel.load(pjoin(self.get_config("DATA_path"), peaks_filename))
            DIRs = niiPEAKS.get_data().astype(np.float32)
            nDIR = np.floor(DIRs.shape[3] / 3)
            print "\t* peaks dim = %d x %d x %d x %d" % DIRs.shape[:4]
            if DIRs.shape[:3] != self.niiMASK_img.shape[:3]:
                raise ValueError("PEAKS geometry does not match with DWI data")

        # setup other output files
        MAPs = np.zeros(
            [
                self.get_config("dim")[0],
                self.get_config("dim")[1],
                self.get_config("dim")[2],
                len(self.model.maps_name),
            ],
            dtype=np.float32,
        )

        if self.get_config("doComputeNRMSE"):
            NRMSE = np.zeros(
                [self.get_config("dim")[0], self.get_config("dim")[1], self.get_config("dim")[2]], dtype=np.float32
            )

        if self.get_config("doSaveCorrectedDWI"):
            DWI_corrected = np.zeros(self.niiDWI.shape, dtype=np.float32)

        # fit the model to the data
        # =========================
        t = time.time()
        progress = ProgressBar(n=totVoxels, prefix="   ", erase=True)
        for iz in xrange(self.niiMASK_img.shape[2]):
            for iy in xrange(self.niiMASK_img.shape[1]):
                for ix in xrange(self.niiMASK_img.shape[0]):
                    if self.niiMASK_img[ix, iy, iz] == 0:
                        continue

                    # prepare the signal
                    y = self.niiDWI_img[ix, iy, iz, :].astype(np.float64)
                    y[y < 0] = 0  # [NOTE] this should not happen!

                    if self.scheme.b0_count > 0:
                        b0 = np.mean(y[self.scheme.b0_idx])

                    if self.get_config("doNormalizeSignal") and self.scheme.b0_count > 0:
                        if b0 > 1e-3:
                            y = y / b0

                    # fitting directions
                    if peaks_filename is None:
                        dirs = DTI.fit(y).directions[0]
                    else:
                        dirs = DIRs[ix, iy, iz, :]

                    # dispatch to the right handler for each model
                    MAPs[ix, iy, iz, :], DIRs[ix, iy, iz, :], x, A = self.model.fit(
                        y, dirs.reshape(-1, 3), self.KERNELS, self.get_config("solver_params")
                    )

                    # compute fitting error
                    if self.get_config("doComputeNRMSE"):
                        y_est = np.dot(A, x)
                        den = np.sum(y ** 2)
                        NRMSE[ix, iy, iz] = np.sqrt(np.sum((y - y_est) ** 2) / den) if den > 1e-16 else 0

                    if self.get_config("doSaveCorrectedDWI"):

                        if self.model.name == "Free-Water":
                            n_iso = len(self.model.d_isos)
                            x[-1 * n_iso :] = 0

                            # print(y, x, b0, A.shape)
#.........這裏部分代碼省略.........

    def fit( self ) :
        """Fit the model to the data iterating over all voxels (in the mask) one after the other.
        Call the appropriate fit() method of the actual model used.
        """
        if self.niiDWI is None :
            raise RuntimeError( 'Data not loaded; call "load_data()" first.' )
        if self.model is None :
            raise RuntimeError( 'Model not set; call "set_model()" first.' )
        if self.KERNELS is None :
            raise RuntimeError( 'Response functions not generated; call "generate_kernels()" and "load_kernels()" first.' )
        if self.KERNELS['model'] != self.model.id :
            raise RuntimeError( 'Response functions were not created with the same model.' )

        self.set_config('fit_time', None)
        totVoxels = np.count_nonzero(self.niiMASK_img)
        print '\n-> Fitting "%s" model to %d voxels:' % ( self.model.name, totVoxels )

        # setup fitting directions
        peaks_filename = self.get_config('peaks_filename')
        if peaks_filename is None :
            DIRs = np.zeros( [self.get_config('dim')[0], self.get_config('dim')[1], self.get_config('dim')[2], 3], dtype=np.float32 )
            nDIR = 1
            gtab = gradient_table( self.scheme.b, self.scheme.raw[:,:3] )
            DTI = dti.TensorModel( gtab )
        else :
            niiPEAKS = nibabel.load( pjoin( self.get_config('DATA_path'), peaks_filename) )
            DIRs = niiPEAKS.get_data().astype(np.float32)
            nDIR = np.floor( DIRs.shape[3]/3 )
            print '\t* peaks dim = %d x %d x %d x %d' % DIRs.shape[:4]
            if DIRs.shape[:3] != self.niiMASK_img.shape[:3] :
                raise ValueError( 'PEAKS geometry does not match with DWI data' )

        # setup other output files
        MAPs = np.zeros( [self.get_config('dim')[0], self.get_config('dim')[1], self.get_config('dim')[2], len(self.model.maps_name)], dtype=np.float32 )
        if self.get_config('doComputeNRMSE') :
            NRMSE = np.zeros( [self.get_config('dim')[0], self.get_config('dim')[1], self.get_config('dim')[2]], dtype=np.float32 )

        # fit the model to the data
        # =========================
        t = time.time()
        progress = ProgressBar( n=totVoxels, prefix="   ", erase=True )
        for iz in xrange(self.niiMASK_img.shape[2]) :
            for iy in xrange(self.niiMASK_img.shape[1]) :
                for ix in xrange(self.niiMASK_img.shape[0]) :
                    if self.niiMASK_img[ix,iy,iz]==0 :
                        continue

                    # prepare the signal
                    y = self.niiDWI_img[ix,iy,iz,:].astype(np.float64)
                    y[ y < 0 ] = 0 # [NOTE] this should not happen!

                    if self.get_config('doNormalizeSignal') and self.scheme.b0_count > 0 :
                        b0 = np.mean( y[self.scheme.b0_idx] )
                        if b0 > 1e-3 :
                            y = y / b0

                    # fitting directions
                    if peaks_filename is None :
                        dirs = DTI.fit( y ).directions[0]
                    else :
                        dirs = DIRs[ix,iy,iz,:]

                    # dispatch to the right handler for each model
                    MAPs[ix,iy,iz,:], DIRs[ix,iy,iz,:], x, A = self.model.fit( y, dirs.reshape(-1,3), self.KERNELS, self.get_config('solver_params') )

                    # compute fitting error
                    if self.get_config('doComputeNRMSE') :
                        y_est = np.dot( A, x )
                        den = np.sum(y**2)
                        NRMSE[ix,iy,iz] = np.sqrt( np.sum((y-y_est)**2) / den ) if den > 1e-16 else 0

                    progress.update()

        self.set_config('fit_time', time.time()-t)
        print '   [ %s ]' % ( time.strftime("%Hh %Mm %Ss", time.gmtime(self.get_config('fit_time')) ) )

        # store results
        self.RESULTS = {}
        self.RESULTS['DIRs']  = DIRs
        self.RESULTS['MAPs']  = MAPs
        if self.get_config('doComputeNRMSE') :
            self.RESULTS['NRMSE'] = NRMSE