import pysam

def calcul_depth_seq(bamfile, seq, window_size=100):

    samfile = pysam.AlignmentFile(bamfile, "rb")
    seq_length = samfile.lengths[samfile.references.index(seq)] #Permet d'obtenir la longueur de la sequence

    for pos_start in range(0, seq_length, window_size): #Itère sur la sequence par "fenêtre" en fonction de la taille
        pos_end = min(pos_start + window_size, seq_length) #pos_start + window_size donne la position de début de la fenêtre, puis grâce à min() et seq_length je peux définir position de fin. le min me permet d'être certain que pos_end < seq_lenght. C-a-d que si la fenêtre dépasse la taille de la sequence, j'aurais quand même pos_end
        window_depths = [] #liste pour stocker la profondeur de chaque fenêtre

        for pileupcolumn in samfile.pileup(seq, start=pos_start, stop=pos_end): #itère sur chaque colone de pileup
            window_depths.append(pileupcolumn.n)#ajoute le nb de lecture couvrant cette position à la liste window_depths

        if window_depths:#calcul la profondeur moyenne pour chaque sequence
            avg_depth = sum(window_depths) / len(window_depths)
            print(f"Profondeur moyenne pour la fenêtre {seq}:{pos_start}-{pos_end}: {avg_depth}")
           
    samfile.close()
    
def calcul_depth_all(bamfile, window_size=100):
    samfile=pysam.AlignmentFile(bamfile, "rb")
    
    for seq in samfile.references:
        calcul_depth_seq(bamfile, seq, window_size)
    samfile.close()

bamfile_path = "/home/theo/dev/git/cnvCallerGPU/align.bam"
#sequence_of_interest = "chr21"
calcul_depth_all(bamfile_path, window_size=100)