package main;
import java.io.BufferedReader;
import java.io.File;
import java.io.FileNotFoundException;
import java.io.FileReader;
import java.io.IOException;
import java.util.NoSuchElementException;
import java.util.Vector;

import javax.swing.UIManager;
import javax.swing.UnsupportedLookAndFeelException;

import org.biojava.bio.*;
import org.biojava.bio.seq.*;
import org.biojavax.bio.seq.RichSequence;
import org.biojavax.bio.seq.RichSequence.IOTools;
import org.biojavax.bio.seq.RichSequenceIterator;
import org.jcvi.jillion.core.datastore.DataStoreProviderHint;
import org.jcvi.jillion.core.qual.PhredQuality;
import org.jcvi.jillion.core.qual.QualitySequence;
import org.jcvi.jillion.core.qual.QualitySequenceBuilder;
import org.jcvi.jillion.core.util.iter.StreamingIterator;
import org.jcvi.jillion.trace.chromat.Chromatogram;
import org.jcvi.jillion.trace.chromat.ChromatogramFactory;
import org.jcvi.jillion.trace.fastq.FastqDataStore;
import org.jcvi.jillion.trace.fastq.FastqFileDataStoreBuilder;
import org.jcvi.jillion.trace.fastq.FastqQualityCodec;
import org.jcvi.jillion.trace.fastq.FastqRecord;

import utils.CompareSequence;
import utils.KMERLocation;
import utils.Subject;

public class TestSingleSequenceQuality {

	public static void main(String[] args) {
		RichSequence seq = null;
		System.out.println(CompareSequence.getOneLineHeader());
		BufferedReader is = null;
		
		//FILL IN
		File singleFileF = new File("Z:\\Datasets - NGS, UV_TMP, MMP\\Targeted Sequencing\\Hartwig\\20210402_PacBio\\raw\\J-Schimmel\\CCS.G09--G09.bam.fq"); 
		FastqDataStore datastoreF = null;
		StreamingIterator<FastqRecord> itF = null;
		try {
			datastoreF = new FastqFileDataStoreBuilder(singleFileF).qualityCodec(FastqQualityCodec.SANGER).hint(DataStoreProviderHint.ITERATION_ONLY).build();
			itF = datastoreF.iterator();
		} catch (IOException e2) {
			// TODO Auto-generated catch block
			e2.printStackTrace();
		}
		
		
		FastqRecord F = null;
		while(itF.hasNext()) {
			F = itF.next();
			if(F.getId().contains("m64072_210430_164756/524304/ccs")){
				break;
			}
		}
		System.out.println(F.getId());
		QualitySequence qs = F.getQualitySequence();
		for(long i=0;i<qs.getLength();i++) {
			if(qs.get(i).getErrorProbability()>0.05) {
				System.out.println(i+"\t"+qs.get(i)+"\t"+qs.get(i).getErrorProbability());
			}
		}
		
		
		String refFile = "Z:\\Datasets - NGS, UV_TMP, MMP\\Targeted Sequencing\\Hartwig\\20210402_PacBio\\ref\\HPRT-exon3_PacBio_o1540-o1541.txt"; 
		//String seqS = " ";
		//String qual = "C%IIHJI:VNBJAWI@VAW@VDJ-WWVFHIII@WJKE9G6HHHKJJJHAVJHAV@CIHI,>IHKHHIII<RJIFH)7=8.SFAW@V'DJHEHHHJ?GEIJJIHC7+VWFIII63=7G<VLH?VJH9WJHAVI(&GIHAVJI@VJF<VJG@U6D.OH?3UKHAVF4CHHJ&=B)&38H@VH?6<@3G3GLJH.JGAVAWHAVJG0@2A75D6D.PJ>L4JIHH6F/TJI@VJH:VJI0VFE8IJ7C-FJ=E;@VKHHKJHI>H=HJ8II&@TJJIB+E9GHF4HIH2;E?VI7H:VJ$>GJII&VUL3H:TWHG>JJ:WVHKJHHJJH4IIIHJ9W?UI@V@7H1JH@WJHB0C8I*NOPOJ<50=%K>(;+5FIHHJKJHHJJIH?H=HJJIIACIHJH@T6IHH/S?VJIH;@W@WJJI5WK9GJJJI>TJ7HHHKG'UGGHHID2QAWE4CI,6PAW:VVA0@H7HHH=E;H?E>HKJGH2BGII'4>G/:IIHH6F?<B*<D7FADI??HA@778-9=GJHH25EG2./0.'%/G:/E:7667FGGIIGH;I5H7H?,87?79HGHJIJ'7CGGOEMGHQKDIRSOLBH+KSQK'QUQN?UQQMVQQSVQQRURQB=F>4)A=:FKJODOQR?OEH?3DMORSPPQJMJ?8NF8QUPRSUPRSVOORV4LPV7FFNF67K;OMTQRSRPLRU@STVSST=fIcEeQ5e\\<QSRD@R-dbSUTOV3^QKTQDPGdOMEaS@FRRFfUTVPCJEH@HfR6ASUE[9HcUOJ->H09><RS&`SBQ@;$AQ&K?DB?849OCRGVQ:+LT=@8)9=.:2=)&CE`SERA'JG9d>L8BBHEM:6DL3P*G9@%P8J7;?WQ7:D5_F`NRU@feQTPSEcRLcD`OEFC_?f_6OGcRQTU:ffeTTTSTHb;.E:VPEN:`KHCeP>BP-`(_Q>>ZB=E<b^0USURK8=?GVOTBQ4XTEEUQUOGRTUE[BccQSUO:)Q=LS:_RT:aNFSUIcUD`QUSWNCQQLD`7^_V?;cNP*OS7e^eUPUQN@%URQVPFfHcR>D5YM3M6ZI3ZAFaSSV3\\ELPSTUUSTR.UPU=RHfTUTTJeT*_=J@29J>OOOOL&:9DMMOMNLLDNK=(4<::D0TQS=U3ebAO5QLOOKMSGeUTHfURUQVJ2PGdRRVT%`TPMKfIc+_QT6RTYcTUP9D5XTTJdITQ@a]-QDMT7eb:XKIJY7e_I^SO;eR&941QM;eQ'ZIVV:aUTU127=BS=VT>AFZQVTVSJeSVGeU:aQGeU>LSQQUSTVPJeS>5aUASQFOQSDZ:DSHfS2ffeQTIeTSL_GfVRTTTJe8efbUR*aa@de@eV>KDeQTUQSKbS7SPTAM0^G.acdbTTSAAfUEHB6<[PUPPV,PKeRTS9eb@VOV7\\LUKVP2ceeP?UTVTVPU76D(@M>MSTVQT%a]^RMQKHdJeR+&DJMZcccXQPRL9RRVTSTFU:DMSJcTS6eeQ3fg^DRPAUQTQUGZSRH:<OSEb+T3bfQ8eeecR:b`SO<=]Q:edSCPTSUTT*^6YbJWQURST@fH(bBPS/fffRIcTFaN-RH@F9]LPW@;;TSUPJe3eefeeSWJQTSUTSO>?.$ZbbbaUD;GGL<<&]JeRUP:SSSOQ9ObcP=VPUGd@SRDQPRRWTSSWRSHcId@AbcPVQDe7eb=TTHfU?bUTV:fe[Q8TOTUTTIeST?GeVSUGf+P,99aOVVS8YB%?XK:cNQBRVSW-CCLUTSPPSTSHb>R<D[TTATR'c`N3L?MgV,e\\SNFW?LQRDKRS4deV;dTNSTLUOQREQ:bSTTTQQUP:edQST<`8PSUO1b&_`YU=IcFfH_8dff]GI>DPK8RSH0ZQP-YRPTUT<R;K=6R;ac0$09^TUTSIc@@9dEY&XN3BRUQSTRC?VVHaU%bcbUM(bbQPJ=BbeV6eee[GPIH51eJ=<LK&;OHeKRSSV?TO2:]PPJaSO=MBHSUO&````TARIOFW$NF;_dePS$4BEJeRCb+@Z;feURK;dV/e_URUQB7]]PUO;dQRV@ffTQA=XeQTT7J_SFLT>XUQJcQFfVPCCd@LJeSVRQ=OO5XRM(bbVPQS?NAQSL+_TUQQRQ:U?OUGf>R9aVVC\\BTSSWPEUR9dVO@cRPUSALJe@R=FX7?<cFeQSQ8fKTRU6bPUTSU5gfQ=RTT9eZ(_NR'Sb^W=cZRPDK8ddRUQP>QTTTSSK?QT1BJOKURU3R:ePAKIQIPJVT@Y:eebDHXRFWTURTSUTHfTGfS>RSM817CI9A_JeFfPL,SQUPOJSVPVS9O?TP9NQRU8JSU:cF`TVN6Z9eeW=;2IbV:ffeGcR5jSSTDAM06IT>P`Q<SQEO:1X\\,*$;J;:2Zff?(RHHfTL;Q2P0DAJ%:a<M:>QVPQVT6gVVO2PFeQ7gVVPLTAgfSHQNFQTN<RVTLBVFfPS6b;YB;<?^bUB.TJe?L'@NQSSJeDfRNHbV?,a[RKQ&4ZeeGe2d_A,Q,7C?J7@KMBeeRJeTRSHU.^JHQ>SIcQUPQUN<HeTTHeG_L8JC\\F9dRV8`NPP4QM<M*F_`RNS@H.UHfSVT?G_5`D7:*OUSVF_AcaPVTTS88;V@efUPFfTJeBS+<N=`TCX@bU4dffaT>ffQVPQSTS>Bb]RU'H[R;aO)_TQ7RacZSNA@IL5V5OEPRTGb>fdD3[feJeQSSAKCJDQVQGX0c`QJcQVT+ccUQWQHbP?\\TC2C?dbUPTTU<ee\\?VVPQTSSV:ff_*D7A7`Dd8e[@HcRLbTVPPWVO+5ObVPNVNLEeS9beVQD?c\\7&GQ5bbbb\\4VG^4R9OWVUSU:fZJeTVRN)__WQPKfRT9aH`SVP9_\\LTJe?NTTR?ff@efQ,8FWSCW:dP5KBT8dEWFPS;cSTFP?P2W\\HcPMG4deTQENDI7LOFd<?3_VHNCVTSTGeB?MD\\J'`^'%\\PVP:UU%bbbZRTTTSGUF`Eb*cccVPSU<V5_DTI4\\E;bQSS9eRVVFTOTHfR<JbOQHTBVE`P=99UG]PFQTTT8`PRQQTJeTUSVR6`V?ffUTHcSFcGD+UL>cBI)\\JRETSWT26I2>T^QWTJe5<TR4[UTVFePJ4.7TU?Je;deQ=GP;0ZffSNJeS@SSA@Q>IdTUHbTRIaQUPOD?deKe.X3)Zaa^0RUQ6_G8_VQIcFe?^H]VTR<bFUOISSIbO7edYG;dTUTBeeSTVTGbITQSQT0PAF>9G<MTSBbdR5ddPALUPVT3eeaRTTVS<?3J]:cSR3deb]SU>URD3eY*aa9JCceTVP-IUcOGfSS@<-CCLIPTT/B51SOOSUT2efeRTSWJPU([`aPFfCQHIdTIeTHfUTSRUM&J]AacQURUHcKfVOSQP?CbFV@H_SS='BX<fcUMS6gN:bVPFdQ:TIeS%SU%JPQ6Q0fcOKP/`T3LYM2[=A5,1LeeVG^PMF@JMVRD_TUJcKWP7eeTRECDaQ=5aP5E'`bbJ-Y@;bOR5_:OMORSBcRHeUVPTUVJeSU@eeIeTS%```^UVK=RJd:b[I3T=PK6M8cOUGVEQ6PHbPCTRK=$HKDeK`3dX5,baGYKeVPA]54DXPF@JJXQQ2NOQRP9_P@STKSEKSJKBTOB0AMM3a;K1f_BTSRUTTUV5cVSRQFdSSTSHcSTVPQSQ>TTTSTSD\\Q1OBDZS>;6]FdR-]PPKFAdSORT?\\O?=LH;KSQ?TSUPGfWC@1SAccUS@SPBGO7US:ccWPTSRBRJcQQRSO@MBRVFTDR>Fd9bUP:$RSEeFb>G`=R/P``RQGd1.1UbTUPVV;ER?OUEaPKeHbUIcQSFc>NALUPR6fFWS9?^eQCeTT3eQN=cedR7L]RSGR='B]T/Wc[TT9eeecR8cdRS+`\\TTQ(+TNR*WaOV3eeT>TIeGfV:LHR-bbVPQJ3J^=BbeTEQIeTVFfQSQ'eedFWJeTT(bbaKGLaS23`edQT7ebQTRVVPLaPB5DBKCO3ad\\QVR'ab`TJ.GX=T0bbVT&a^PHQ]aaa^RS8FVdcSWT+dd^KH?VP)b\\MAO^cQ2A0WF=.`=TTJcT%12DX^M:dEeTCeeIc:PUdQP;eebP7c>;>LJe7cTN4de]S@GQPBecR@edTWTJVM9_eRVJe?QVR6feTNCeT9aUEIOJ5J0&UL+gggc_UR*\\bQGMT:_ebU(bbU:beeU;eeaC>`S9J3]edHe;=UGdQ1dee]Q9F_?S2+BD<,DD@B2;GKS+^CEWQ4Uce_R:d`SBebK>HTAE?DBeeS6^QPT;dQSSRDC@.3:^_R>7eQ7I]5dfS)cddddbT9aI1cb<aecGfSIeTV<eeePS=G`C8eeeaGfS$`^X_\\L7/3D3aaQ1UTDfPJISOVPIF=JeDaPPWTSSFNTPISJ;bN:;\\F+*;1;9_V?5`bQQT)ZZZZYU4`RTUR?CCTBeeSUQ4eebR4`aQWTFGSTINQOPSOCD;GSPTQTPQIRNL584<>RQTQTQUPV7eZBcT:bDcSNS6[0ddd^QQVPCecRP>-_dG>>.DM\\ZQBecP.fffSUVPM3aR.``VPSJcSTVNQN>.a_NL3*K\\cW@RO?S6fD,I,QbdX%THCAR9ffaU*^_ab]URTP(^_YGdU-d]M=<;AHDL6SaaaR6XR9dXEV)abbTQAe_A2WD;0ZbK:VPUSPGMSVTVR=eedSU6]SQ+V>)_S5eZO=8Z`SN/J;beH?9Z_FfQSJd?8`ZLSQS,cdddddOGfVNSQQVPCLTIeRQ<SCee,bZQ;>,beeeeQ?eb97b]RKP<`QT4``R4OOSUQ.ZcddJ=IeTT6^edRS-LB<J:*C8OAVS,caDA:2F]e`Q)^LFT6L\\BK:K3AJO4]eeRQS.aaa[QTSARVFfPWTVBeTTTJeR3cebRHc+bcbWLP?T?1K;H8H.7]YTVP<d\\K6cdQTUQPBSJeG_'\\]YFSD+WP6`SUVOR;MPVVP,ccc]P@HbTPI(^RTS0cS]PKRHGO=fQPTI_7cJeSKeQEf:d_TQ8F-\\\\MBQ2JMQ@KOO6^cN5ceVAeeT8dU=5^bR&)HdDgT9KeEQB:deTUTT;e_SO4cedR9dWL,LbdT7WJeSVQLB+YSR<e[NT=rss^Fs]OZPL\\Qua\\QtBcLj_[\\Oo[Ip[>cdPs\\F>nteQ0_]R=]o],MCI+pqqqk]]L\\`^Pu_`_\\Qu\\\\NT^:ntsMaQ>l[_^BG^\\5q`]Qt\\[%leVZ(pp[.9ZkX)kpqqncJD\\Q]_\\Qu\\[^Duj[8X_Y;O83iqq\\FV\\AtZ^^\\UA7Qt]GGsrCtt_Mj0@p]Pn]3`W0NJN]B@gTQq_5_N(;H]^El]IZLBU]IW]Po^NrP/ss_\\]L\\D'dprrrr_]Qs\\\\S\\PbNYACncI@m_UV8cqt\\L_FZ\\\\_\\_\\9uut^VBncUZEa[HQg<MGCY^\\^\\]_]_`XO&booTQs[\\OsX,YQsSo^2rtth_`_\\.tte_a\\9qiBRh_\\:~J6SE$.Rqtg]_1ttdPsRu_\\WJ_RoK0nsd]_4tp_g^Pu]Iq]E\\9rr^L^^\\X3bB[]^Rs^[KMtJ;X]<RqZ9rpl])oqqmIdJ^EoJpIAst\\[HurLXX6][Js][?:lSm\\_OUKNu\\Rs\\XAo`\\-iWaU\\QCn]_]ZORF;qnCq\\6svk__^_\\V`ZG[[__Z_@saAtut[4se^[^]]?R]\\Gtp[*T^\\K\\\\?N\\`W]=lX\\_LjALpT(eppS_]\\`_[]<k]_Qv[Qt^^[CC,]ck^]FY&pqe,bTPRHttO`Qs^AuiQsQsZBr[Qs2CHO[^][NUQs\\[;lPrY]I`SqEGps^<tu_[\\:pLZ_^J][_\\^`Lq_[^QsKRYU@u\\a\\\\]M^QK\\\\Qu_^Dsu_a_\\Qs_QWDY\\[_`\\Cuu^IMnLY\\J]DtX@i\\0rOBp^^__PtXX\\D\\YJEaQ'J]K]KY&eoX^Pu_a[Os^[V]\\+robhAqt\\4V*Q,Zqtttts]\\@i^[_7nNe]\\[KXIECRZI]\\_RtQs8vnZ6Wdtm^\\\\`Qs]\\\\Cus^[OoI^JW^[H_^]_`YKhF]I_\\]PqH;^W_T@hXOc\\__\\^^\\_Pu^^FrpZ7[^JXC^WC_]]:qtmI^^JJT]IX__\\\\_Iut[[_`Y_KNAr>K1ZROc\\]2@m`]_^[XJ]_^\\_MQ^W_8MbPt\\^]4mlG)245&]RFJ@q\\_Os^QsPtPdJ[[_GT<UVNs^Bup^O^_\\_\\Duu\\[_5sPB_UJDIRJMLHOR?ELIQBSUP6MMPUTSQ5*5Q\\YSIKL=GEPNTPE\\[[K^]JX\\MML]^1hsaQo2rk\\E]_Qh\\_\\_UGPG>qqY[_`\\7tts_`\\7tut^X_Os[_^Ns_\\\\\\Qsa_^\\_@[BKIg2M8?UF\\]M?;rs\\La'Rcpqqn]_C\\[AuLXQ^BLErrJqX8eW4V_QsNsEFo";
		String seqS = F.getNucleotideSequence().toString();
		QualitySequence qual = F.getQualitySequence();
		String left = "AAGACTTGCTCGAGATGTCA";  
		String right = "TGAAGGAGATGGGAGGCCAT"; 
		//String left = null;
		//String right = null;
		
		File ref = new File(refFile);
		try {
			//is = new BufferedReader(new FileReader("Z:\\Joost\\Files\\Manuscripts\\Schimmel_etal_2016\\Robin\\HPRT-FASTA-CR2.txt"));
			//is = new BufferedReader(new FileReader("Z:\\Joost\\Files\\Manuscripts\\Schimmel_etal_2016\\Robin\\HPRT-FASTA-CR1.txt"));
			//is = new BufferedReader(new FileReader("Z:\\Evelyne\\JavaSoftware\\XF1290.fa.txt"));
			//is = new BufferedReader(new FileReader("Z:\\Evelyne\\DNA\\Revertants sequencing\\XF1423_extended.fa.txt"));
			//is = new BufferedReader(new FileReader("Z:\\Tim\\G23 insertion\\XF1426.fa.txt"));
			//is = new BufferedReader(new FileReader("C:\\Users\\rvanschendel\\Documents\\Project_Primase\\polq-1_reversion\\XF1280_whole_unc-22_100bp_zone_for_polq.fa"));
			//is = new BufferedReader(new FileReader("E:\\Project_Hartwig\\HPRT-FASTA-CR1.txt"));
			is = new BufferedReader(new FileReader(refFile));
			
		} catch (FileNotFoundException e1) {
			// TODO Auto-generated catch block
			e1.printStackTrace();
			System.exit(0);
		}
		//get a SequenceDB of all sequences in the file
		RichSequenceIterator si = IOTools.readFastaDNA(is, null);
		RichSequence hprtSeq = null;
		while(si.hasNext()){
			try {
				hprtSeq = si.nextRichSequence();
				//IOTools.writeFasta(System.out, hprtSeq, null);
			} catch (NoSuchElementException e) {
				// TODO Auto-generated catch block
				e.printStackTrace();
			} catch (BioException e) {
				// TODO Auto-generated catch block
				e.printStackTrace();
			} catch (Exception e) {
				// TODO Auto-generated catch block
				e.printStackTrace();
			}
		}
			Subject subjectObject = new Subject(hprtSeq, left, right);
			subjectObject.setPacBio(true);
			CompareSequence kmerWith = new CompareSequence(subjectObject,seqS,qual, null, true, "");
			//CompareSequence(RichSequence subject, String query, QualitySequence quals, String left, String right, String pamSite, String dir, boolean checkReverse, String queryName, KMERLocation kmerl) {
			//s.setAdditionalSearchString(additional);
			kmerWith.setAndDetermineCorrectRange(0.05);
			kmerWith.maskSequenceToHighQualityRemove();
			kmerWith.determineFlankPositions(false);
			String outputKMER = kmerWith.toStringOneLine("dummy");
			System.out.println(outputKMER);
			
			/*
			
			kmerless1 = new CompareSequence(hprtSeq, seq.seqString(),chromo.getQualitySequence(), left,right, null, null, true, seq.getName(), null);
			mkerless1.setAndDetermineCorrectRange(0.05);
			//s.setAdditionalSearchString(additional);
			//s.maskSequenceToHighQuality(left, right);
			kmerless1.maskSequenceToHighQualityRemove();
			kmerless1.determineFlankPositions();
			outputKMERless = kmerless1.toStringOneLine();
			System.out.println(outputKMERless);
			
			kmerWith = new CompareSequence(hprtSeq, seq.seqString(),chromo.getQualitySequence(), left,right, null, null, true, seq.getName(), kmerl);
			kmerWith.setAndDetermineCorrectRange(0.05);
			//s.setAdditionalSearchString(additional);
			//s.maskSequenceToHighQuality(left, right);
			kmerWith.maskSequenceToHighQualityRemove();
			kmerWith.determineFlankPositions();
			outputKMER = kmerWith.toStringOneLine();
			System.out.println(outputKMER);
			
			if(kmerWith.getRemarks().length()==0 && kmerless1.getRemarks().length()==0 && !outputKMERless.equals(outputKMER)) {
				System.err.println("KMER AND KMERLESS ARE NOT THE SAME (QUALTIY MASKED)");
				System.err.println(outputKMERless);
				System.err.println(outputKMER);
			}
			*/
	}
}
